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ABSTRACT 

We have modelled the process of reionization of the intergalactic medium (IGM) by 
photoionization by galaxies, in order to learn what galaxy formation in the frame¬ 
work of the CDM model predicts for the epoch of reionization. We use a sophisticated 
semi-analytic model of galaxy formation to track the formation of these galaxies, their 
influence on the IGM, and the back-reaction of the state of the IGM on further galaxy 
formation. Our study represents a much more complete and physically consistent mod¬ 
elling of reionization than has been conducted in the past. In particular, compared to 
previous work by ourselves and others, our new calculations contain signihcant im¬ 
provements in the modelling of the effects of reionization of the IGM on the collapse of 
baryons into dark matter halos (this is now computed self-consistently from the prop¬ 
erties of model galaxies), and in the model for the cooling and condensation of gas 
within halos (our new model includes photoheating from a self-consistently computed 
ionizing background and also includes cooling due to molecular hydrogen). We hnd 
that reionization can be achieved by z ~ 10-20 in a AGDM cosmological model with 
as ~ 0.9. However, a cosmological model with a running spectral index is only able to 
achieve reionization before z ~ 9, and thus be consistent with an optical depth of 0.1, 
if very extreme assumptions are made about the physics of feedback at high redshifts. 
We also consider the specific galaxy formation model recently discussed by Baugh et 
ah, which includes a top-heavy IMF in starbursts, and hnd that it is able to reionize 
the Universe by z w 12. The previous results assume that all of the ionizing photons 
produced by stars in galaxies are able to escape and ionize the IGM. If this is not the 
case, then the redshift of reionization could be substantially reduced. We hnd that 
extended periods of partial reionization and double reionizations can occur in models 
in which the hrst stars formed via cooling by H 2 molecules, are very massive, and 
in which the escape fraction of ionizing photons ^10-30%. Such models do not fully 
reionize until z ~ 6-7, but predict an electron scattering optical depth as large as 0.15. 

Models with lower as = 0.7-0.8 as suggested by the recent Wilkinson Microwave 

Anisotropy Probe (WMAP) three year data have reduced redshifts of reionization, 
but can be consistent with the lower optical depth also suggested by the WMAP three 
year data. 


1 INTRODUCTION 

While gas in the intergalactic medium (IGM) is highly ion¬ 
ized in the present day Universe, it is known to have been 
mostly neutral in the past. Evidence for this comes from 
observations of the Gunn-Peterson effect in the spectra of 
high redshift quasars (Becker et al. 2001; Djorgovski et al. 
2001; Fan et al. 2006) and from the value of the optical 
depth to scattering measured from the cosmic microwave 
background (GMB; Spergel et al. 2003). Thus, there was a 
period in the history of the Universe when it transitioned 
from being neutral to ionized—the epoch of reionization. 

Measurements of the GMB polarization made during 
the first year of operation of the Wilkinson Microwave 


Anisotropy Probe (WMAP) satellite^ have been used to es¬ 
timate a value r — 0.17 ± 0.04 for the electron scattering 
optical depth to reionization (Kogut et al. 2003), which, 
in the simple case that the IGM reionized instantaneously, 
implies a reionization redshift Zreion = 17 ± 5 (Spergel et 
al. 2003). This reionization redshift is significantly higher 
than was previously thought on the basis of the observed 
Gunn-Peterson effect and on the basis of theoretical mod- 


^ While this paper was being completed, results based upon the 
first three years of data from WMAP became available (Spergel et 
al. 2006). These show a much lower optical depth to reionization. 
We will comment upon these results in §5. 
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elling (see, for example, Haiman & Loeb 1997; Gnedin & 
Ostriker 1997; Baugh et al. 1998; Valageas & Silk 1999; Br- 
uscoli et al. 2000; Chiu & Ostriker 2000; Ciardi et al. 2000; 
Benson et al. 2001). Although the WMAP measurement 
could be explained in scenarios in which reionization hap¬ 
pened twice, or happened only partially for some time (Cen 
2003; Wyithe & Loeb 2003; Naselsky & Chaing 2004),^ the 
simple interpretation (that the Universe was fully reionized 
at 2 = 17 and remained so ever since) implies the pres¬ 
ence of a significant population of ionizing sources at high 
redshift. Coupled with the hints that the effective spectral 
index of the matter power spectrum is changing with scale 
(i.e. a “running spectral index”), this scenario becomes hard 
to explain in the standard cold dark matter (CDM) scenario 
of structure and galaxy formation.® 

This tension has been explored recently by several au¬ 
thors (Haiman & Holder 2003; Somerville, Bullock & Livio 
2003; Yoshida et al. 2003b; Ciardi, Ferrara & White 2003; 
Onken & Miralda-Escude 2004), who generally reach the 
conclusion that Zreion = 17 is incompatible with the running 
spectral index model, and hard to achieve even in a power- 
law model. One possible escape from this apparent problem 
is to reionize the Universe in some other way, i.e. without re¬ 
lying on ionizing photons from stars in high redshift galaxies. 
One such mechanism which has received attention recently 
is a contribution to reionization from particle decays during 
the cosmic dark ages (Chen & Kamionkowski 2004; Hansen 
& Haiman 2004; Kasuya, Kawasaki & Sugiyama 2004; 
Kasuya & Kawasaki 2004), another is a contribution from 
“mini-QSOs” (see, for example, Ricotti & Ostriker 2004). 

In this paper, we undertake detailed calculations of 
the galaxy formation process at high redshift to investi¬ 
gate whether galaxies alone can ionize the Universe suffi¬ 
ciently early to explain the WMAP measurement of the op¬ 
tical depth. For this purpose we employ the GALFORM semi- 
analytic model of galaxy formation (Cole et al. 2000; Benson 
et al. 2002a). The advantage of our model over previous at¬ 
tempts to address this question is that it incorporates much 
more of the actual physics of galaxy formation (e.g. it follows 
the cooling of gas within haloes and incorporates a model of 
star formation and feedback), even if this physics is treated 
in a simplified way. The model used in this paper is signifi¬ 
cantly improved over that employed in our earlier study of 
reionization in Benson et al. (2001). The effects of reioniza¬ 
tion on subsequent galaxy formation (due to the heating of 
the IGM and the subsequent rise in the Jeans mass, and due 
to photoheating of gas in haloes) are now self-consistently 
computed within our model. Furthermore, we now include a 
detailed treatment of cooling due to molecular hydrogen in 
our calculations. Thus, we are able to provide a more reliable 
answer as to what the CDM model of structure formation 
predicts for the epoch of reionization. 

The remainder of this paper is arranged as follows. In 

^ The physical plausibility of these scenarios has been questioned 
by Furlanetto & Loeb (2005). We present evidence, in §3.3 that 
such double or partial reionization may in fact occur under plau¬ 
sible physical conditions. 

® The running spectral index reduces the amplitude of fluctua¬ 
tions on small, galactic scales relative to that on large scales. As 
a result there are many fewer haloes present at high redshift with 
mass sufficient to form galaxies. 


§2 we outline the model used to follow the reionization of 
the Universe. In §3 we examine the reionization histories 
predicted by our models and assess the importance of the 
various physical ingredients of the models. In §4 we present 
results for the reionization epoch and the imprint of reion¬ 
ization on the cosmic microwave background. Finally, in §5 
we assess the implications of our results. Several Appendices 
describe further details of the modelling, and present addi¬ 
tional results. 

2 MODELLING 

2.1 Galaxy Formation Modelling 

We employ the GALFORM semi-analytic model of galaxy for¬ 
mation throughout this work to predict the properties of 
the galaxy population as a function of redshift. Most impor¬ 
tantly for the current work is the ability of this model to 
compute the star formation rates (and associated luminosi¬ 
ties in ionizing photons) needed to predict the reionization 
history of the Universe. A full description of the GALFORM 
model is beyond the scope of this paper. We instead refer 
the reader to previous papers which document the model in 
more detail (Cole et al. 2000; Benson et al. 2002a). Here, 
we give a very brief description of key components of the 
model and provide somewhat more detail about two ingre¬ 
dients which are critical to the current work: the cooling of 
hot gas in haloes (see §2.2) and star formation (see §2.3). 
It should be noted that the major innovation with respect 
to Paper I is our improved model of gas cooling, which now 
includes cooling due to H 2 and which is affected by feedback 
from the reionization process. The implications of this new 
modelling are discussed in §2.2. 

The basic structure of the GALFORM calculation can be 
summarized as follows: 

• The hierarchical formation of dark matter haloes is fol¬ 
lowed using the extended Press-Schechter theory (Bond et 
al. 1991; Bower 1991). This allows the construction of 
merger trees which form the backbone of the galaxy for¬ 
mation calculation. 

• A galaxy can form in each branch of the merger tree if 
the hot gas in the halo is able to cool and condense to the 
halo centre (see §2.2) where it forms a rotationally supported 
disk. 

• This condensed gas is assumed to begin to form stars 
at a rate governed by simple, phenomenological rules (see 
§2.3). As a consequence of this star formation, supernovae 
inject energy into the galaxy causing some of the condensed 
gas to be re-ejected from the galaxy. This feedback loop is 
crucial to producing a population of galaxies which matches 
the observed low-redshift Universe (Benson et al. 2003a). 
Star formation also produces metals, enriching the gas in the 
galaxy. When this enriched gas is ejected from the galaxy by 
supernova explosions, the hot atmosphere of gas pervading 
the halo also becomes enriched in metals. 

• When two branches of the merger tree join, the most 
massive galaxy from the branches becomes the central 
galaxy of the new branch. Any other galaxies become satel¬ 
lite galaxies, free to orbit within the gravitational potential 
of the new dark matter halo. Satellite galaxies are no longer 
able to gain gas from the surrounding hot atmosphere. 
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• Satellite galaxies experience dynamical friction as they 
orbit, which eventually forces them to merge into the galaxy 
residing at the halo centre. These mergers are assumed to 
transform stellar disks into stellar spheroids if the masses 
of the two galaxies are comparable. Thus, in this model, 
mergers are the drivers of morphological evolution. Mergers 
of this type may also trigger bursts of star formation (see 
§2.3). 

2.2 Cooling of Gas in Haloes 

When a dark matter halo reaches virial equilibrium, it is 
assumed to have accreted into itself some mass of gas from 
the surrounding IGM. For high mass haloes, the mass of gas 
accreted is equal to f2b/(f^o — fib) times the dark matter 
mass of the halo (where flo and fib are the mean densities 
of total and baryonic mass in the Universe in units of the 
critical density). For lower mass haloes, the mass of gas ac¬ 
creted is reduced, due to the finite gas pressure in the IGM, 
in accordance with the filtering mass formalism of Gnedin 
(2000). The filtering mass, Mp, is discussed in detail in Ap¬ 
pendix A3. The mass of gas accreted by a halo of mass Mtot 
is then given by 

j./r _ _ ( fib/flp) Mtot _ ...s 

S- [l + (2l/3_l)MF(2)/Mtot]3- 

The gas is assumed to be heated by shocks as it ac¬ 
cretes into the halo, reaching a uniform temperature equal 
to the halo virial temperature. The gas is then assumed to 
be distributed through the halo with a density profile 

f's-W « (2) 

where c is a scale radius which we set equal to rs/3, where 
Ts is the scale radius of the dark matter density profile (as¬ 
sumed to have the NFW form). 

To compute the rate of cooling of the gas, we con¬ 
struct a table of net radiative cooling rates as a function 
of temperature, density, metallicity and redshift. (The red- 
shift dependence arises through the inclusion of photoheat¬ 
ing from a time-varying ionizing background which is com¬ 
puted self-consistently from the emission from the model 
galaxies.) This net cooling/heating rate includes contribu¬ 
tions from atomic processes (e.g. recombinations, collisional 
excitations etc.), Compton cooling (Peebles 1968), cooling 
due to molecular hydrogen (see §A1) and photoheating by 
the ionizing background (see §2.4). Given these cooling rates 
and the assumed gas density profile and temperature for 
each halo, we then calculate, at each time step of our calcu¬ 
lation, the cooling radius, defined as the radius in the halo at 
which the local cooling time equals the time since the halo 
formed. (A halo is defined to form anew whenever its mass 
is doubled by mergers.) Any gas within this cooling radius is 
assumed to accrete onto the galaxy forming at the centre of 
the halo, provided that the time for the gas to free-fall to the 
halo centre is less than the age of the halo. Occasionally, for 
example when the ionizing background is rapidly rising, the 
cooling radius at a given timestep may be smaller than that 
at the previous timestep. In such cases no gas is allowed to 
cool until the cooling radius once again begins to increase. 

The inclusion of cooling due to H 2 is crucial to the 
present work as, at high redshifts, this is the only cooling 


channel available to gas in the majority of halos. Without 
H 2 cooling included, we risk seriously underestimating the 
epoch of reionization. The IGM pressure and the ionizing 
background act as sources of negative feedback on star for¬ 
mation (i.e. star formation leads to a rise in the IGM pres¬ 
sure and ionizing background, which in turn suppress later 
star formation). By including these effects self-consistently 
in our model we are therefore able to check whether this 
negative feedback is sufficient to significantly quench star 
formation. It should be noted that if we neglect H 2 cooling 
and set the filtering mass and ionizing background to zero 
our revised cooling model produces results identical to those 
of Paper I (i.e. it simply follows the cooling prescription of 
Cole et al. (2000) with the addition of Compton cooling). 
We will examine the effects of including H 2 cooling in our 
calculations in §3.1. 

The effects of including the filtering mass and ionizing 
background in our cooling model have been explored in de¬ 
tail by Benson et al. (2002a). Benson et al. (2002a) find 
that it is the filtering mass which has the greatest effect 
on galaxy formation at low redshifts, resulting in a sup¬ 
pression of star formation after the epochs of Hi and Hel 
reionization. For the Cole et al. (2000) model, Benson et 
al. (2002a) found this suppression to be quite small. How¬ 
ever, for a model with weak feedback from supernovae the 
suppression effect of the filtering mass was significant (re¬ 
ducing the global star formation rate by a factor of two). 
This occurs because, in models with weak SNe feedback, 
star formation can otherwise progress efficiently in low mass 
halos, and therefore such haloes make a significant contri¬ 
bution to the global star formation rate. Once the filtering 
mass rises star formation is shut off in these low mass halos, 
thereby greatly suppressing the global star formation rate. 
The filtering mass also acts to reduce the number of faint 
galaxies which can form, as discussed in detail by Benson et 
al. (2002b). Note that the filtering mass begins to rise as the 
process of reionization begins, but does not reach its maxi¬ 
mum value until significantly after the epoch of reionization. 
As such, we may expect it to have only minimal impact on 
the epoch of reionization, a point which we will explore in 
§3.2. 


2.3 Star Formation and Supernova Feedback 

Once gas has cooled inside a dark matter halo it is as¬ 
sumed to collapse, conserving its angular momentum, until 
it reaches a radius where it becomes rotationally supported 
against further collapse. At this point, the gas is assumed to 
form an exponential disk. The gas in this disk then begins 
to form stars at a rate ^|) given by 

!/> = Mcoid/r*, (3) 

where Mcoid is the total mass of gas available in the galaxy 
and r* is a characteristic timescale for star formation. Fol¬ 
lowing Cole et al. (2000) we adopt the parameterization 

r* = e* Vdisk(V'disc/200 km/s)“*, (4) 

where e* and a* are dimensionless parameters, Vdisc is the 
circular velocity of the disk at its half-mass radius and rdisc 
is the dynamical time of the disk at that radius. 
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We assume that supernova explosions in the galaxy 
cause gas to be ejected from the galaxy at a rate Meject = 

where 

/?= (5) 

where Ohot is a dimensionless parameter and Vhot is a pa¬ 
rameter with dimensions of speed. 

Star formation can also occur in a burst mode as the 
result of a violent merger between two galaxies. When such 
a merger occurs, all of the gas in the two merging galaxies 
undergoes star formation on a very short timescale. We use 
the same rules as for quiescent star formation in disks, but 
using the dynamical properties of the newly formed spheroid 
instead. 

We compute the chemical enrichment of each galaxy us¬ 
ing the instantaneous recycling approximation, with yield of 
heavy elements p and recycled fraction R {as defined by Cole 
et al. (2000)—we use the same values for these parameters 
as did Cole et al. (2000) unless otherwise noted), includ¬ 
ing the exchange of metals between the galaxy and the halo 
gas. Once the star formation history and associated chem¬ 
ical enrichment history of a galaxy have been determined, 
we compute the SED of the galaxy at each timestep using a 
stellar population synthesis model based on the Padova stel¬ 
lar evolution tracks (Granato et al. 2000). This allows us to 
compute the ionizing photon luminosity for each galaxy. 


2.4 Evolution of the IGM and Ionizing 
Background 

We model the process of reionization and its feedback on the 
IGM and galaxy formation using the methods of Benson et 
al. (2001; hereafter, Paper I) and Benson et al. (2002a). 
Briefly, we compute the volume averaged emissivity in ion¬ 
izing photons at each redshift by summing the ionizing lu¬ 
minosity of all galaxies in our calculation. This emissivity is 
used to evolve a distribution of gas elements drawn from an 
appropriate PDF (Benson et al. 2001). We track the ion¬ 
ization (H and He) and thermal state of these gas elements 
as the Universe evolves. From the resulting thermal history 
of the IGM we are able to compute the filtering mass de¬ 
scribing the effects of the IGM pressure on the collapse of 
baryons into haloes. This allows us to estimate the mass of 
baryons accreted into each dark matter halo forming at any 
redshift. We also have, at each redshift, the spectrum of the 
ionizing background, which allows us to compute the rate of 
photoionization heating experienced by gas in dark matter 
haloes. These two processes (reduced accretion due to the 
IGM pressure and photoheating of halo gas), which act to 
suppress galaxy formation, are incorporated into a second 
iteration of our galaxy formation calculation. This iterative 
process is repeated until convergence is reached (i.e. the in¬ 
put ionizing background and filtering mass are consistent 
with those resulting from the galaxy formation calculation) 
as in Benson et al. (2001). Schematically, the iterative pro¬ 
cedure works as follows: 

(i) Initially assume that there is no ionizing background 
(and associated photoheating) and that the filtering mass is 
always zero. 

(ii) Run the GALFORM galaxy formation calculation. 


(iii) Compute the evolution of the ionizing background 
based on the GALFORM calculation. 

(iv) Compute associated photoheating rates for gas of 
varying densities, temperatures and metallicities. 

(v) Compute the evolution of the filtering mass based on 
the GALFORM calculation. 

(vi) Using these photoheating rates and filtering mass, go 
back to step (ii) and repeat. 

Iteration stops once the star formation history produced in 
GALFORM calculations on successive iterations is sufficiently 
well converged. 

Benson et al. (2001) and Benson et al. (2002a) con¬ 
sidered cooling due to both atomic processes and Comp¬ 
ton cooling due to the scattering of cosmic microwave back¬ 
ground photons from hot electrons. Here we employ an im¬ 
proved cooling model which is described in detail in Ap¬ 
pendix A. To explore the high redshifts and low mass scales 
of interest here we must also consider the effects of cooling 
due to molecular hydrogen (see §A1). We must also recon¬ 
sider Compton cooling in the non-equilibrium regime (i.e. 
before free electrons in the gas have had sufficient time to 
recombine to reach their abundance in ionization equilib¬ 
rium at the density of the halo gas) which occurs in some of 
the lower mass haloes that we consider (see §A2). 

2.5 Escape Fractions 

A very important parameter in our models is the fraction, 
/esc, of ionizing photons produced in a galaxy which es¬ 
cape into the IGM. For most of our models (and unless 
stated otherwise), we assume /esc = 1. The only excep¬ 
tions are the models H 2 WeakFB and H 2 WeakFB-VMS (to 
be defined in §2.7.2), for which we also present some re¬ 
sults for lower escape fractions (for the purpose of producing 
models with “double reionization”—see §3.3). The assump¬ 
tion that /esc = 1 is obviously an extreme one (although 
such an assumption has been used in other works studying 
reionization, e.g. Haiman & Holder 2003), since it ignores 
the absorption of ionizing photons by gas and dust in the 
galaxy, but it should result in the earliest possible reioniza¬ 
tion if other parameters are held fixed. We note, however, 
that observational evidence (and simple physical models for 
/esc) imply a much lower escape fraction. Therefore, in Ap¬ 
pendix B we present a simple, physically motivated calcula¬ 
tion of /esc which is in reasonable agreement with observa¬ 
tional data, and show the resulting reionization histories. 

2.6 Calculation of Filling Factor 

To compute the filling factor of Hll we follow the growth 
of ionized regions around each halo in the galaxy formation 
model using the methodology of Paper I. Specifically, we 
compute the total emission rate of ionizing photons from 
stars in each dark matter halo, S{t), and compute the num¬ 
ber of hydrogen atoms which will be ionized by this radiation 
by solving 

= S{t) - a^^'a'^/ciumpUHiVH, (6) 

where A/j is the total number of hydrogen atoms ionized, 

is the case B recombination coefficient for hydrogen (for 
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a temperature of lO'^K), a is the cosmic expansion factor, 
/clump is the clumping factor for ionized gas, and nn is the 
mean number density of hydrogen atoms in the Universe. We 
compute the clumping factor following the method described 
in Paper I for , but with the cosmological parameters 

appropriate to the current work. In this model, the main 
contribution to the clumping factor comes from gas residing 
in haloes with virial temperatures around lO^K, which have 
deep enough potential wells to retain gas, but which are 
cool enough to be not completely collisionally ionized. The 
specihc values for the clumping factors used in this work will 
be discussed in §2.7.1. 

Note that we do not assume any particular geometry 
for the ionized region around each source. Ionized regions 
from two or more sources may overlap (particularly just be¬ 
fore the epoch of reionization or if the sources are strongly 
clustered). This does not affect our results, as we are simply 
counting the total number of hydrogen atoms which have 
been ionized. 

The total number of hydrogen atoms, Nn.tot, which will 
have been ionized in some representative volume of the Uni¬ 
verse, V, is then found by simply summing over the contri¬ 
butions from galaxies in haloes of all masses 

A^h ,tot — E Nn.iniU, (7) 



M/h 'Mg 

Figure 1. The fractional root mean square fluctuation of the 
density averaged in spheres, a{M), as a function of mass M, for 
each of our cosmological models: WMAP RSI (red line), WMAP 
PL (blue line) and Baugh et al. (2005; green line). 


where Nii,i is the number of hydrogen atoms which have 
been ionized by galaxies in halo i, of which there are Ui such 
haloes per unit volume. The sum is taken over all haloes. 
The mean ionized fraction in volume V is found by simply 
dividing by the total number of hydrogen atoms in the vol¬ 
ume, nnU where nn is the mean density of hydrogen atoms 
in the Universe. If we assume a uniform IGM, we can identify 
this mean ionized fraction with the filling factor of ionized 
regions, 

Tfiii = — Nu.irii. (8) 

riH ^ 
i 

We consider the Universe to have been fully reionized once 
Tflii = 1. Note that with these dehnitions Tfui can exceed 
unity. In such cases, the Universe is fully ionized and ionizing 
photons will begin to build up an ionizing background rather 
than be absorbed by neutral hydrogen atoms. 

The IGM evolution model of Benson et al. (2002a) 
gives an alternative means of following the reionization of 
the Universe. In that model, the mean ionization state of 
the IGM is determined by computing the ionization and 
recombination rates in the presence of an evolving ionizing 
background, allowing for a distribution of IGM densities. 
While these two methods for following reionization are quite 
different, we find that they predict very similar redshifts 
of reionization, with the method of Benson et al. (2002a) 
typically predicting ^reion to occur roughly Az « 1 earlier. 
Throughout this work we will show results based on the 
approach described in this section. 

2.7 Cosmological and Galaxy Formation Models 

We will explore three different cosmological models and sev¬ 
eral different models of galaxy formation in this work. The 
details of these models are outlined below. 


2.7.1 Cosmological Parameters 

For the cosmological parameters we adopt the values de¬ 
termined from the hrst-year WMAP data set (Spergel et 
al. 2003). In particular, we consider cosmological models 
with the best ht parameters for a power-law ACDM model 
with constraints from WMAP alone (Table 7 of Spergel et 
al. 2003) and a running spectral index AGDM model with 
constraints from WMAP, the 2dF Galaxy Redshift Survey 
and Lyman-a forest observations (Table 8 of Spergel et al. 
2003). We refer to these cosmological models as WMAP PL 
and WMAP RSI respectively, and their parameters are given 
in Table 1. We also compute reionization histories in the cos¬ 
mological model corresponding to the specihc galaxy forma¬ 
tion model of Baugh et al. (2005), discussed further below. 

We include these three cosmologies because they have 
complementary strengths and weaknesses. The cosmologies 
based on WMAP results give good flts to the available CMB 
data (by deflnition). However, we do not have a set of galaxy 
formation model parameters in these cosmologies which re¬ 
sult in galaxy populations with properties which match local 
observations. The Baugh et al. (2005) model, on the other 
hand, was not designed to match the GMB data (although 
its cosmological parameters are not too far from those mea¬ 
sured by WMAP), but does reproduce the properties of 
galaxies in the local and high-redshift Universe. 

We show in Figure 1, for all three cosmological mod¬ 
els, the fractional fluctuation of the density measured in 
spheres, a{M), for the linear density held extrapolated to 
z = 0. The quantity a{M) determines the abundances of 
dark matter haloes at different redshifts. It is clearly seen 
that the WMAP RSI model has much smaller fluctuations 
on small scales. 

The IGM clumping factor which we use depends only 
on the cosmological model, and is shown as a function of 
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Table 1. Cosmological parameter sets used in this work. Values are taken from (Spergel et al. 2003; WMAP) and (MacTavish et al. 
2005; Boomerang), ris is the scalar spectral index at A: = 0.05Mpc“^, while dris/dlnfc is the running of the spectral index as defined by 
eqn. (5) of Spergel et al. (2003). 


Name 

Uq 

Ao 

Ub 

i7o/km s 1 Mpc ^ 

0-8 

ns (at k = 0.05Mpc 

dns/d In k 

WMAP PL 

0.270 

0.730 

0.0463 

72.0 

0.90 

0.99 

0.000 

WMAP RSI 

0.268 

0.732 

0.0444 

71.0 

0.83 

0.93 

-0.031 

Baugh et al. ’05 

0.300 

0.700 

0.0400 

70.0 

0.93 

1.00 

0.000 

Boomerang ’03 

0.300 

0.700 

0.0469 

69.6 

0.85 

0.95 

0.000 



Figure 2. The clumping factor, /dumpi El'S a function of redshift 
for the three cosmological models considered in this work. Also 
shown, for reference, are the clumping factors corresponding to 
the cosmological parameters of Paper I (see Appendix Cl) and 
the “Boomerang ’03” cosmology discussed in §5. 

redshift in Fig. 2. For the WMAP PL cosmology it is very 
similar to that found in Paper I (/dump « 18 at 2 = 10, drop¬ 
ping to approximately 2.5 at z = 20 for example), while for 
the WMAP RSI cosmology it is much lower (/dump ~ 1-5 at 
a = 10, dropping to close to 1 at z = 20) due to the reduction 
of small-scale fluctuations in the latter cosmology. Inclusion 
of the clumping factor makes almost no difference to our 
results for the WMAP RSI cosmology, as the clumping fac¬ 
tor is still small by the time reionization occurs. However, 
in the WMAP PL cosmology, clumping can delay reioniza¬ 
tion by Az = 2-4 (the earlier reionization occurs, the less 
effect clumping has). In the Baugh et al. ’05 cosmology, the 
clumping factor is comparable to that of the WMAP PL 
cosmology. 

2.7.S Galaxy Formation Parameters 

To model the formation of galaxies in the Universe we em¬ 
ploy the semi-analytic model of galaxy formation GALFORM. 
A full description of this model can be found in Cole et 
al. (2000) and Benson et al. (2002a). For our “Standard 
model” in the WMAP RSI and WMAP PL cosmologies, we 
adopt parameters which are very similar to those of Cole et 


al. (2000). (We ignore the cosmological parameters speci¬ 
fied by Cole et al. (2000) and instead use those described in 
§2.7.1.) However, we use the improved merging calculation 
of Benson et al. (2002a), which results in somewhat fewer 
mergers on average than the simpler model used by Cole et 
al. (2000). Due to the change in cosmological parameters 
and in the merger rates, we adjust the values of a few of 
the parameters in the galaxy formation model, increasing 
the star formation efficiency e* from 0.005 to 0.0067 and 
reducing the mass-to-light ratio parameter T from 1.38 to 
1 . 00 . 

We regard the resulting model as being a plausible 
model for galaxy formation in the high-redshift Universe, 
since it has been shown to reproduce well the properties 
of galaxies in the low-redshift universe, and we simply as¬ 
sume that the same physical processes apply to galaxy for¬ 
mation at high redshift as at low redshift. Our model is 
therefore motivated by what we currently know about the 
local Universe—the extrapolation to high redshifts is a sim¬ 
ple assumption as our current understanding of the physics 
that governs the high redshift Universe does not allow us to 
devise more accurate physical prescriptions for early times. 
Note, however, that the Cole et al. (2000) model assumed 
Db = 0.02. The higher baryon densities in the cosmologies 
considered here will result in this model producing too many 
galaxies of all luminosities at z = 0. We do not attempt to 
correct for this deficiency for two reasons: 1) we are inter¬ 
ested here in determining how early reionization can occur— 
overproduction of galaxies will push reionization to higher 
redshifts and so give us a conservative answer to the ques¬ 
tion of whether a particular model can reionize sufficiently 
early to match the WMAP results; and 2) while the Cole et 
al. (2000) model works well at a = 0 it does not match ob¬ 
servational data for high-redshift galaxies and, in any case, 
we are here extrapolating our model to « > 20 where ob¬ 
servational constraints on galaxy formation are practically 
non-existent. It is for these reasons that we consider our 
model to be plausible but not definitive. 

While there are numerous variations of the model pa¬ 
rameters which it would be interesting to explore, we limit 
ourselves here to examining a handful which might be ex¬ 
pected to have a particularly large impact on the epoch of 
reionization. These various models are summarized in Ta¬ 
ble 2. 

• Standard: For the WMAP RSI and WMAP PL cos¬ 
mologies, this model is very similar to that described in Cole 
et al. (2000)—the differences from Cole et al. (2000) are the 
parameter changes described in §2.7.2 and the inclusion of 
Compton cooling, which is necessary for studying the high 
redshifts of interest here. It assumes that all stars form with 
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Figure 3. The mass of a halo with virial velocity corresponding to 
the feedback parameter, Vhot = 200km/s, is shown as a function 
of redshift for the WMAP PL cosmology (results for the WMAP 
RSI and Baugh et al. ’05 cosmologies are practically identical). 


a Kennicutt (1983) initial mass function (IMF), with slope 
X = 0.4 for 0.1 < m < IMq and x = 1.5 for 1 < m < 
12OM0. (A Salpeter IMF has x = 1.35.) For this IMF, the 
total number of ionizing photons (with hv >13.6 eV) pro¬ 
duced per Mq of stars formed is A^uyc/Af* = 1.0 x 10®®, 
for a metallicity Z = O.1.^0. Note that, for this model 
(and all others unless stated otherwise) we adopt a value 
of Vhot = 200km/s for the feedback parameter. This value, 
chosen to reproduce the low abundance of faint galaxies in 
the z = 0 Universe, results in quite strong feedback. For the 
Baugh et al. (2005) cosmology, the Standard model corre¬ 
sponds to the galaxy formation parameters given by Baugh 
et al. (2005). Note that this includes a higher value for 
Vhot = 300km/s. Figure 3 shows the halo mass with virial 
velocity equal to Vhot = 200km/s as a function of redshift. 
At high redshifts, the vast majority of haloes able to form 
will be strongly affected by feedback for Vhot ^ 200km/s. 

• MetFree: In this model, we compute stellar luminosi¬ 
ties and spectra assuming that all stars are be metal free 
(while retaining the chemical enrichment of the interstel¬ 
lar and intergalactic gas, as this affects gas cooling rates), 
thereby enhancing the emission of ionizing photons from 
galaxies. In this case, we have Ahyc/Af* = 2.0 x 10®®. 

• TopHeavy: We assume an extreme IMF which has the 
Salpeter slope x = 1.35, but only covers the mass range 10 < 
m < IOOM 0 , i.e. it is truncated for masses below IOM 0 . As 
massive stars produce the majority of ionizing photons, this 
model is expected to produce earlier reionization. For this 
IMF, we have N-^yc/Mt = 7.1 x 10®®. 

• WeakFB: This model has very weak feedback from su¬ 
pernova explosions, such that even the lowest mass galaxies 
which form in the model will not experience signihcant mass 
loss due to feedback. In GALFORM, the strength of feedback 
is controlled by two parameters Vhot and Ohot, as dehned 
in eqn.(5). We fix Ohot = 2 throughout this work, and use 


Vhot = 200(300)km/s for the WMAP (Baugh et al. ’05) cos¬ 
mologies except for models with weak feedback, for which 
we assume Vhot = Ikm/s (see Table 2). The lower value of 
Vhot = Ikm/s corresponds to an efficiency of supernova en¬ 
ergy injection and of mass ejection which is lower by a factor 
4 X 10"^ than for the standard case Vhot = 200km/s, so it is 
fairly extreme. This weak feedback will greatly enhance the 
efficiency of star formation in low-mass haloes, and hence 
will also increase the production of ionizing photons at high 
redshift. 

• H 2 : This model is identical to the Standard model 
except that it incorporates cooling due to molecular 
hydrogen—this should allow galaxies to form earlier, thereby 
shifting the epoch of reionization to higher redshifts. 

• H 2 WeakFB: This model combines a weak feedback 
with molecular hydrogen cooling. 

• H 2 WeakFB-VMS: This model investigates the possi¬ 
bility that stars forming via H 2 cooling (“Population III”) 
may be very massive stars (VMS), with m ~ 100 — IOOOM 0 
(Abel, Bryan & Norman 2000; Abel, Bryan & Norman 
2002; Bromm, Coppi & Larson 2002). If all of the stars 
are very massive, then the total ionizing luminosity will 
be greatly enhanced relative to stars formed with a normal 
IMF. Therefore, in this model, we increase the UV luminos¬ 
ity of any stars forming in a halo with virial temperature be¬ 
low lO^K by a factor of 20 relative to the Kennicutt IMF at 
all metallicities.® We therefore have Ni^yRMt = 2.0 x 10®^. 

Baugh et al. (2005) recently presented a GALFORM 
model employing a significantly different set of galaxy for¬ 
mation parameters from those used by Cole et al. (2000). 
We will use the Baugh et al. (2005) parameters to de¬ 
scribe our “Standard model” in the Baugh et al. (2005) 
cosmology. The cosmological parameters of this model are 
listed in Table 1. The Baugh et al. (2005) model was con¬ 
strained to produce a good match to observations of sub-mm 
galaxies and Lyman-break galaxies at high redshifts, as well 
as the properties of galaxies at low redshifts. In order to 
match the numbers of sub-mm galaxies, Baugh et al. as¬ 
sumed that stars which form in merger-induced bursts have 
a top-heavy IMF (with slope a: = 0 covering the mass range 


^ Our choice of a UV enhancement factor of 20 is motivated by 
the results of Bromm, Kudritzki & Loeb (2001), who find an 
enhancement factor of 10—20. (Their enhancement is relative to a 
Salpeter IMF. We use the Kennicutt (1983) IMF here, but the to¬ 
tal ionizing photon production per unit mass of star formation for 
this IMF differs from that of the Salpeter IMF by approximately 
10% for sub-Solar metallicities.) We have taken the value at the 
upper end of their range to maximize the effects. We also note 
that theoretical considerations suggest that very massive stars 
should form when the metallicity of the ISM is below a critical 
value, and not merely when gas cooling occurs through the H 2 
channel. If we model VMS formation in this way (i.e. we enhance 
the UV emissivity only for stars with a metallicity below, for ex¬ 
ample, 10“'* Z 0 (Schneider et al. 2002; Bromm & Loeb 2004) 
we find little enhancement in the net UV emission. This occurs 
because the duration of the sub-critical metallicity phase of star 
formation is very short. However, our current model does not ac¬ 
count for the different heavy element yield of these Population 
III stars—a low yield could dramatically extend the duration of 
the sub-critical metallicity phase of star formation. We defer a 
detailed examination of this possibility to a future work. 
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0.1 < m < I 2 OM 0 ), while stars which form quiescently in 
disks have a Kennicutt (1983) IMF (as in our Standard 
model). In the Baugh et al. model, the fraction of stars 
formed in the burst mode increases with redshift, with the 
burst mode dominating the total star formation density at 
2 > 3 (see Fig. 1 in Baugh et al. ). In subsequent papers, 
Le Delliou et al. (2005a,b) have investigated the properties 
of Lya-emitting galaxies at high redshift, and Nagashima et 
al. (2005a,b) have investigated the chemical enrichment of 
elliptical galaxies and intracluster gas in the same model. 
Given the success of the Baugh et al. model in reproducing 
observations of various types of star-forming galaxies at high 
redshift, as well as observations of chemical abundances at 
low redshift, it is of interest to see what it predicts for the 
reionization of the IGM. We have therefore used the galaxy 
formation and cosmological parameters of the Baugh et al. 
(2005) model to compute the reionization history of the Uni¬ 
verse® . 

For stars formed in bursts in the Baugh et al. model, 
the number of ionizing photons produced per Mq of stars 
formed A^Lyc/Af* = 1.7 x 10®^ for a metallicity Z — O.IZq, 
i.e. nearly 20 times larger than for a Kennicutt IMF. The 
only model variants we consider for the Baugh et al. cos¬ 
mology at Standard, H 2 and H 2 WeakFB®. 


3 RESULTS 

Below, we present results from our calculations of reioniza¬ 
tion. We begin with a summary of the key results, before 
concentrating on the importance of H 2 cooling and pho¬ 
toionization feedback in setting the epoch of reionization. 
A comparison of the predicted epochs of reionization to the 
observational determinations and an exploration of the sig¬ 
nal imprinted in the cosmic microwave background (GMB) 
by the patchy reionization process are given in §4. 

Some further results, including comparisons to Paper I 
and tests of the convergence of our models are given in Ap¬ 
pendix G. Appendix B briefly examines the effects on the 
reionization history of using physically motivated calcula¬ 
tions of the escape fraction of ionizing photons. 

Figure 4 shows key results from all of our calculations. 
Each row corresponds to a different model (as indicated 
along the right-hand edge of the figure) and within each 
panel different colours correspond to different cosmologies: 
WMAP RSI (red), WMAP PL (blue) and Baugh et al. ’05 
(green). The three columns respectively show the average 
star formation rate per unit volume, the ionizing luminosity 
density and the Hll filling factor (as indicated by labels at 
the top of the figure) all as functions of redshift. It is clear 
from this figure that the suppression of small mass haloes in 

® Note that Baugh et al. (2005) included a simple model of 
photoionization feedback in their calculation, in which they sup¬ 
pressed the cooling of gas in dark matter haloes with virial veloci¬ 
ties below 60km/s for z < 6. As we are here modelling reionization 
and photoionization feedback in much greater detail, we remove 
this simple suppression from our implementation of the Baugh et 
al. (2005) model. 

® Model H 2 WeakFB-VMS is not considered, as the Baugh et al. 
(2005) model already includes very massive stars which are as¬ 
sumed to form in bursts. 



Figure 5. The Hll filling factor as a function of redshift for models 
with (dot-short dashed lines) and without (solid lines) H 2 cooling. 
Colours indicate the cosmological model as indicated in the figure. 
Results are shown to higher redshifts and lower filling factors than 
in Fig. 4. 

the WMAP RSI cosmology leads to a severe reduction in the 
star formation rate, ionizing luminosity and Hll filling fac¬ 
tor compared to our other cosmological models. This effect 
is larger than that of any of the other variations we consider 
in this work. As such, it may be possible to use a measured 
reionization epoch to place strong constraints on the cos¬ 
mological model. It is also clear that reionization redshifts 
high enough to explain the WMAP measurement of the op¬ 
tical depth can be attained in the WMAP PL and Baugh 
et al. ’05 cosmologies, albeit only in rather extreme models. 
These points will be discussed further in the remainder of 
this section and in §4. 

3.1 Does H 2 Cooling Help? 

We can examine whether the inclusion of H 2 cooling chan¬ 
nels in our calculations actually aids in pushing the epoch 
of reionization to higher redshifts. Figure 5 shows the Hll 
filling factor as a function of redshift for models with and 
without cooling due to H 2 . It can clearly be seen that the 
addition of H 2 cooling channels has almost no effect on the 
epoch of reionization. At earlier redshifts, H 2 cooling does 
result in an enhancement in the filling factor (e.g. by a factor 
of 10^ at 2 = 35 in the WMAP PL cosmology). However, by 
the time Ffui is approaching unity, the ionizing emissivity is 
dominated by emission from haloes which are able to cool 
through atomic processes, and so the addition of H 2 cooling 
channels does not help to increase 2reion significantly. 

To understand the effects of cooling through molecu¬ 
lar hydrogen, we can explore the cooled baryonic fraction 
in haloes as a function of their mass. The upper panels of 
Fig. 6 show the condensed baryonic fractions of haloes at 
redshifts 5 (left-hand panel) and 10 (right-hand panel) for 
models without the effects of photo-ionization and photo- 
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Table 2. Names and details of the galaxy formation models used in this work. Column 1 describes the model while column 2 gives 
the name by which we refer to the model throughout this work. Column 3 lists the changes in model parameters made relative to the 
standard model. 


Model 

Name 

Variation from standard model 

Standard model 

Standard 

— 

Metal free stellar spectra 

Met Free 

2^* = 0 

Top-heavy IMF 

TopHeavy 

IMF Salpeter above IOM 0 , zero below; R = 0.91; p = 0.14 

Weak feedback 

WeakFB 

Vhot = Ikm/s 

Compton -h H 2 cooling and dynamical heating 

H 2 

H 2 cooling on 

Compton -h H 2 cooling -h weak feedback 

H 2 WeakFB 

H 2 cooling on; Vhot = Ikm/s 

Compton -h H 2 cooling -h weak feedback 
+ VMS for Tvir < lO-^K 

H 2 WeakFB-VMS 

H 2 cooling on; Vhot = Ikm/s; UV flux x20 for stars formed 
via H 2 cooling 


Figure 4. (cont.) The evolution of key physical parameters in our models. The three columns show mean star formation rate density, 
ionizing luminosity and Hll filling factor as indicated at the top of each panel. Each row corresponds to a different galaxy formation model, 
as indicated on the right-hand side of the figure. The three lines in each panel correspond to three different cosmological models: WMAP 
RSI (red), WMAP PL (blue) and Baugh et al. (2005; green—not shown for models MetFree, TopHeavy, WeakFB and H 2 WeakFB-VMS). 
In the Hll filling factor columns, the horizontal dashed line indicates the point of complete reionization, while dashed coloured lines for 
the H 2 WeakFB-VMS model show the contribution of Pop I and II stars to the total filling factor. 


dissociation included. Blue circles indicate the cooled frac¬ 
tion, /cond (i.e. the mass of material which has been able to 
cool and condense into the galactic phase), for a model with 
no molecular hydrogen cooling, while red circles show results 
for a model with molecular hydrogen cooling. We indicate, 
with a vertical dashed magenta line, the halo mass for which 
the virial temperature corresponds to the supernova feed¬ 
back parameter Vhot in our semi-analytic model. Vertical, 
dashed cyan lines indicate the filtering masses for the two 
models. In the lower panels of Fig. 6 are curves which show 
the cooling measure tage/(1 + P)tcoo\, the ratio of the age of 
the Universe at the specified redshift tage to the cooling time 
fcooi (calculated at 200 times the mean density of the Uni¬ 
verse, which is approximately the mean density of haloes), 
and divided by 1 -|- /3, where /3 (given by eqn. 5) is the ratio 
of the mass ejection rate due to supernova feedback to the 
star formation rate in model galaxies. Solid blue lines are for 
the model without molecular hydrogen cooling and without 
Compton cooling, while dashed red lines include molecular 
hydrogen cooling (and Compton cooling). We show curves 
for primordial gas and for gas of half Solar metallicity (the 
upper curves at large halo masses). 

Our cooling measure, the quantity tage/(l + /3)tcooi, is 
interesting since it roughly determines whether a halo has 
had time to turn most of its supply of gas into stars. Roughly 
speaking, a galaxy will accrete (from its surrounding hot at¬ 
mosphere) all available gas in a time tcooi. Assuming haloes 
to have survived for a period roughly equal to tage, the ra¬ 
tio of these two timescales then determines whether or not 
a halo has been able to accrete all available gas. However, 
much of this gas will be ejected once again due to feedback. 
This gas may be subsequently re-accreted. Consequently, a 
galaxy must accrete its hot halo ~ {1 + (3) times before it 
can turn the majority of the available gas into stars.^ 

^ Assuming the cooling/reheating loop to occur n times (and 
assuming cooling and reheating to happen instantaneously), the 
fraction of available gas eventually bound into stars will be 1 — 


When tage/(l + P)tcooi is greater than one (the point 
fage/(l + P)tcooi = 1 is indicate by green horizontal lines 
in Fig. 6), the effective cooling time is short compared to 
the age of the Universe and we do not expect the cooling 
function to have strongly determined the shape of the /cond 
curve. When it is less than unity, the effective cooling time 
is long and so the shape of the /cond curve should have been 
determined strongly by the shape of the cooling curve. 

From Fig. 6 we can see that the inclusion of H 2 cool¬ 
ing makes no difference to the fraction of condensed mass 
in haloes with virial temperatures above about 2 x lO^K. 
H 2 cooling does however, result in a tail of low mass haloes 
(with virial temperatures below lO^K) which are able to 
condense some of their gas, unlike in the case of cooling 
through purely atomic processes. However, the amount of 
gas condensed is minimal, with < 10% of the total available 
gas being condensed (a number which falls rapidly to zero 
in lower mass haloes). This inefficient condensation occurs 
because tage/tcooi(l + /3) 1 (typically 10“®-0.1) for these 

haloes due mostly to the fact that /3 S> 1 (typically, (3 > 100 
for these haloes)—this makes star formation a very ineffi¬ 
cient process in these galaxies, with the majority of the gas 
which was able to condense being re-ejected due to feedback. 

Further plots of the condensed fraction for all cosmolo¬ 
gies and models can be found in Appendix D. 

3.2 Effects of including photoionization and 
photodissociation feedback 

We have run all models both with and without the effects 
of photoionization feedback (i.e. the suppression of accre¬ 
tion due to the IGM pressure and photoheating of gas in 
haloes by an ionizing background) and photodissociation 
feedback (i.e. the destruction of H 2 molecules by the dissoci¬ 
ating background). In Fig. 7 we show the resulting star for- 

(/3/(l + /3))". For n = 1 + /3 this converges, for large /3 to 1 — 1/e Ri 
0.63. 
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Figure 6. Measures of the ability of baryonic material to cool into the galactic phase are shown for 2 = 5 (left-hand panels) and 2 = 10 
(right-hand panels) for models without the effects of photo-ionization and photo-dissociation feedback. Vertical dashed magenta lines 
indicate the halo mass for which the virial temperature corresponds to the feedback parameter Vbot- Vertical, dashed cyan lines indicate 
the filtering masses for the two models. Upper panels: The average fraction of baryonic material in haloes which has been able to cool and 
condense into galaxies, as a function of halo mass in the Standard model in the WMAP RSI cosmology. Blue circles indicate the cooled 
fraction for a model with no molecular hydrogen cooling, while red circles show results for a model with molecular hydrogen cooling. 
Lower panels: Cooling measures for metallicities Z = 0 and Z = 0.5Zq (the upper curves at large halo masses): solid blue lines are 
for the model without molecular hydrogen cooling, while dashed red lines include molecular hydrogen cooling (and Compton cooling). 
The cooling measure is the ratio tage/(l + /3)fcool- Here 0 is the ratio of rate of mass ejection due to supernova feedback to the star 
formation rate in model galaxies. The horizontal green lines indicate the point where the cooling measure equals 1 (i.e. the point above 
which cooling is efficient in the sense described in the text. 


mation rates and Hll filling factors as a function of redshift. 
Photoionization and photodissociation feedback reduce the 
redshift of reionization in general due to the suppression of 
galaxy and star formation. Fig. 7 also illustrates once more 
that the inclusion of H 2 cooling does little to increase the 
redshift of reionization. In fact, it appears that including 
H 2 cooling actually reduces the redshift of reionization in 
the curves plotted here. This is a numerical artifact, arising 
from the fact that the most massive haloes simulated in our 
calculations are of lower mass in models which include H 2 


cooling (see §C2.2). Once this is corrected for, reionization 
redshifts are almost identical whether or not H 2 cooling is 
included (see §4.1). 

Figure 8 shows the same information as Fig. 6 but for 
models which now include the effects of photo-ionization and 
photo-dissociation feedback. By comparing this figure with 
Fig. 6, it can be seen from this figure that the main causes of 
the suppression of galaxy formation are photoheating of hot 
gas in haloes and photodissociation of H 2 . The filtering mass 
is approximately the same in models with and without pho- 
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Figure 7. Star formation histories (blue and cyan lines) and Hll filling factors calculated using ionized spheres (red and magenta lines). 
Cyan and magenta lines show results for models where photoionization and photodissociation feedback are neglected, while blue and 
red lines show results for models including these feedback mechanisms. Results are shown for Standard and H 2 models in WMAP PL, 
WMAP RSI and Baugh et al. ’05 cosmologies as indicated for each column and row. 


toionization and photodissociation feedback included® (it is 
not until after reionization that the filtering mass has the 
chance to respond to the increase in the temperature of the 
gas). At z = 5 we can clearly see that the cooling times 
in haloes with virial temperatures below around 10®K are 


® Note that, even when photoionization feedback is not included 
the filtering mass still rises due to the heating of the IGM by 
ionizing photons emitted by galaxies. All that is left out of the 
calculation is the effect of the filtering mass on subsequent galaxy 
formation. 


greatly increased due to this photoheating and photodisso¬ 
ciation. At 2 = 10 the effect is much weaker, due both to 
the lower background at that redshift and the fact that gas 
in haloes is denser at higher redshifts and so less affected by 
photoheating. 

It can be seen that the inclusion of photoionization and 
photodissociation feedback has only a small effect on the 
redshift of reionization. This is not surprising, as the mech¬ 
anisms through which photoionization feedback works (pho¬ 
toheating of gas in haloes by the ionizing background and an 
increase in the filtering mass) become important only after 
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Figure 8. As Fig. 6 but for models including the effects of photo-ionization and photo-dissociation feedback. 


reionization. Figure 9 shows the evolution of the Hll filling 
factor®, filtering mass and ionizing background for the stan- 


® The “spikes” seen in the Hll filling factor are an artifact of the 
semi-analytic model. As in Paper I, we compute the reionization 
process by drawing a set of halo masses from the halo mass func¬ 
tion at some redshift z and constructing merger trees for each 
halo back to very high redshifts. The galaxy formation calcula¬ 
tion then proceeds in each tree. Our merger tree algorithm is 
not completely accurate in the sense that, when used to produce 
merger trees that span large redshift ranges, it does not exactly 
reproduce the expected distribution of progenitor halo masses. 
This problem becomes more severe as the redshift range spanned 
increases. Consequently, we repeat our calculations for several 
starting redshifts Zi (ensuring that at each Zi we have the correct 
distribution of halo masses). The reionization histories shown are 
constructed by splicing together results from each calculation, us¬ 
ing quantities between Zi and from each calculation. At each 
Zi we find discontinuities in the predicted properties due to this 
splicing process. 


dard model in the WMAP RSI cosmology with and without 
the effects of photoionization feedback included. The filter¬ 
ing mass begins to rise shortly before reionization occurs, 
when photoheating first dominates over cooling of gas due 
to the expansion of the Universe, but does not achieve any 
significant value until after the Universe is already reion¬ 
ized. Similarly, the ionizing background cannot build up to 
any significant level (and thereby efficiently photoheat gas 
in haloes) until the Universe has already been ionized and 
thereby become transparent to ionizing radiation. Benson 
et al. (2002a) found that it was the filtering mass which 
caused the greatest suppression of galaxy formation, seem¬ 
ingly in contradiction with the above result. Benson et al. 
(2002a) were, however, considering the effects of photoion¬ 
ization feedback on the galaxy population at a = 0, by which 
time the filtering mass had had time to fully respond to the 
change in the IGM temperature during reionization. When 
considering the effects of photoionization feedback on the 
epoch of reionization itself, as we do here, the filtering mass 
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Figure 9. Indicators of the strength of photoionization feedback 
mechanisms as a function of redshift for the Standard model in the 
WMAP RSI cosmology. Results are shown for this single model 
and cosmology only as they are typical of all of the models we 
consider. The solid red (blue) line in the upper panel shows the 
filtering mass (scale on left-hand axis) when photoionization feed¬ 
back is (is not) included. In the lower panel red (blue) circles show 
the ionizing background intensity (scale on left-hand axis) when 
photoionization feedback is (is not) included. The dashed red line 
in both panels indicates the Hll filling factor (scale on right-hand 
axis) for a model including photoionization feedback. 


is much less important as it does not have time to fully re¬ 
spond to the heating of the IGM until well after reionization 
is over. The ionizing background, however, rises very rapidly 
as reionization occurs, allowing photoheating to have some 
(albeit small) effect on the redshift of reionization. 


3.2.1 Effects of H 2 Dissociation? 

In general, we find that stars formed via H 2 cooling con¬ 
tribute little to the reionization of the Universe, if they form 
with a normal IMF. Therefore, it is not surprising that the 
inclusion or otherwise of H 2 dissociation in our model makes 
no significant difference to the end results on the ionized fill¬ 
ing factor, as shown in Fig. 10. 



Figure 10. The Hll filling factor as a function of redshift for 
model H 2 with (solid line) and without (dashed line) dissociation 
of H 2 molecules included in our calculations. Results are shown 
for the WMAP RSI cosmology only, but are typical of all other 
cosmologies and models. 

3.3 Enhanced Emission from Early Generations of 
Stars 

In model “H 2 WeakFB-VMS” we enhance the UV emission 
from stellar populations which form via H 2 cooling (i.e. those 
which form in dark matter haloes with virial temperatures 
below 10‘^K) by a factor of 20 relative to a standard IMF, to 
account for the possibility that these stars might be super- 
massive “pop III” stars. The left-hand panel of Figure 11 
shows the resulting Hll filling factors vs. redshift for this 
model and for model H 2 WeakFB (which is identical apart 
from having no enhancement of UV emission due to super- 
massive stars) for the WMAP cosmologies. The results are 
plotted for an escape fraction chosen to give a double reion¬ 
ization: 25% (WMAP RSI) and 15% (WMAP PL). This 
results in all models producing full reionization at 2 « 6. 
Note however, that model H 2 WeakFB-VMS displays an ex¬ 
tended period of partial reionization back to higher redshifts 
(from 2 « 8 for the WMAP RSI cosmology to 2 17 

for the WMAP PL cosmology). As a consequence, model 
H 2 WeakFB-VMS has an increased electron scattering opti¬ 
cal depth as shown in the right-hand panel of Figure 11. 
For example in the WMAP PL cosmology we find r = 0.15, 
compared to 0.08 for model H 2 WeakFB. Values of r are 
lower for the other cosmologies, but still come significantly 
closer to the WMAP measurement of r = 0.17 ± 0.04. As 
such, model H 2 WeakFB-VMS is in much better agreement 
with the data obtained from the WMAP satellite (Kogut et 
al. 2003). Enhanced UV emission from stellar populations 
formed via H 2 cooling is therefore able to reconcile the ap¬ 
parent disparity between a large optical depth and evidence 
for reionization at 2 « 6 in the WMAP PL cosmology, by 
inducing a period of partial reionization. 

The dotted line in Fig 11(a) shows the contribution to 
the Hll filling factor from stars forming via H 2 cooling in 
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the H 2 WeakFB-VMS model for the WMAP PL cosmology. 
Before z ~ 20 these stars dominated the reionization of the 
Universe, as the steepness of the dark matter halo mass func¬ 
tion results in there being very few haloes hot enough to cool 
via atomic processes at these redshifts. At lower redshifts the 
abundance of haloes with virial temperature above lO^K in¬ 
creases, while galaxy formation in lower temperature haloes 
becomes inefficient due to increased cooling times and the 
rising filtering mass^°. Therefore, haloes hotter than lO'^K 
rapidly become the dominant contributors to the Hll filling 
factor. 

Obtaining a period of extended partial reionization, or a 
double reionization, requires an interplay of several key fac¬ 
tors. Firstly, the cosmological model must produce enough 
haloes with virial temperatures below lO'^K to permit a sig¬ 
nificant Hll filling factor to be obtained from the early gen¬ 
eration of stars (assumed to be supermassive). The enhance¬ 
ment factor (taken to be 20 in this work) is also important: 
if it is too high then reionization may happen once only, if 
it is too low the early generation of stars will not be able 
to produce a significant Hll filling factor. The decline in the 
contribution from early generation stars as time progresses 
must also be sufficiently rapid (otherwise the filling factor 
will continue to rise, leading to a single reionization event)— 
requiring either that the filtering mass grows rapidly or that 
cooling times for haloes cooler than lO^K increase rapidly 
shortly after the first reionization. Despite this requirement 
that many factors combine to produce a double reionization, 
we find, for quite reasonable values of the UV enhancement 
factor, that double reionization is a natural outcome of our 
H 2 WeakFB-VMS model in all cosmologies considered. An 
extensive survey of parameter space would be required to 
determine how common such double reionizations are, but 
this is beyond the scope of this work. 

While our modelling of this enhanced UV emission from 
early stars is very simplified (e.g. we assume that supernova 
feedback is the same for Very Massive Stars and for normal 
stellar populations), and based on uncertain assumptions 
(e.g. that all haloes which cool by H 2 produce Very Mas¬ 
sive Stars), it is still relatively sophisticated compared to 
previous studies of double reionization (Cen 2003; Wyithe 
& Loeb 2003) to the extent that ours is the only model 
which follows the full hierarchical merging history of each 
forming galaxy and, furthermore, computes the fraction of 
baryons able to cool and the fraction able to turn into stars 
using physically motivated prescriptions rather than simply 
making these fractions parameters of the model. 

4 OBSERVABLE QUANTITIES 
4.1 The Epoch of Reionization 

Table 3 fists the redshifts of reionization in all of the models 
considered (all results shown include the effects of photoion¬ 
ization and photodissociation feedback), along with the cor¬ 
responding electron-scattering optical depths. Observations 
of the cosmic microwave background suggest a large value 

The filtering mass in the H 2 WeakFB-VMS model in the 
WMAP PL cosmology rises above the mass corresponding to a 
lO^K halo at z = 15. 


for Zreion- Kogut et al. (2003) used WMAP measurements 
of the temperature-polarization (TE) cross-power spectrum 
to estimate the optical depth r. They found a best-fit value 
r = 0.17, with a 95% confidence range 0.09 ^ r 4 0.28. For 
the WMAP cosmology, this corresponds to a best-fit value 
of the reionization redshift Zreion = 17, and a 95% confidence 
range 11 < Zreion < 24, if one assumes instantaneous reion¬ 
ization, i.e. if the ionized fraction increased from 0 to 1 at 
2 = 2reion. Spergel et al. (2003) found very similar results 
from their multi-parameter fit to the whole WMAP dataset 
combined with data on galaxy clustering. 

Since our models predict a gradual rather than an in¬ 
stantaneous reionization, we compare our models with the 
WMAP data using the value of r rather than the value of 
Zreion. From Table 3 we see that all of our models based 
on the WMAP PL cosmology and assuming /esc = 1 are 
consistent with the measured r (at 95% confidence), apart 
from the H 2 WeakFB-VMS model, which predicts r = 0.30, 
marginally too high. These models predict 2:reion in the 
range 11-22. However, when /esc is reduced to 15%, the 
H 2 WeakFB-VMS model in this cosmology predicts r = 0.15 
(and Zreion = 17), in good agreement with the WMAP mea¬ 
surement. In contrast, all but one of our models based on 
the WMAP RSI cosmology, and assuming /esc = 1, are ex¬ 
cluded (at 95% confidence) by the WMAP measurement 
of r, due to the predicted optical depths being too low, 
the one exception being the H 2 WeakFB-VMS model, which 
(with r = 0.10) is marginally compatible. The models in 
the WMAP RSI cosmology predict Zreion in the range fi¬ 
ll. Thus, while the WMAP RSI cosmology cannot be com¬ 
pletely excluded based on the WMAP measurement of t, it 
does appear very unlikely. If the escape fraction /esc is in 
fact significantly less than unity, then the predicted reion¬ 
ization epochs will be shifted to considerably lower redshifts 
in all of the models. For the Baugh et al. ’05 cosmology, we 
find that the Standard and H 2 models are marginally con¬ 
sistent with the WMAP constraint, while models with weak 
feedback are able to match that constraint within current 
errors. 

Observations of the spectra of the highest redshift 
quasars at « ~ 6 show evidence of Gunn-Peterson absorp¬ 
tion by neutral hydrogen in the IGM (Becker et al. 2001; 
Djorgovski et al. 2001; White et al. 2003; Fan et al. 2003; 
Fan et al. 200fi), which has been interpreted as evidence 
that the Universe reionized at redshifts only slightly above 
6, since simple theoretical models predict that the transi¬ 
tion from the IGM being almost completely neutral to al¬ 
most completely ionized should have happened over a rather 
narrow redshift range. A low Zreion ~ 7 is in apparent contra¬ 
diction with the results from WMAP. This interpretation of 
the observed quasar spectra has however been questioned by 
some authors (e.g. Songaila 2004). Even if current measure¬ 
ments of the Gunn-Peterson optical depth tgp at 2 ~ 6 are 
correct, they can be explained by very small neutral hydro¬ 
gen fractions in the IGM, e.g. White et al. (2003) measure 
TGP ^5 — 20 for two quasars at 2 « fi.3, which can be pro¬ 
duced with neutral fractions of only /hi ^ 1.5 — fix 10“®. 
Therefore, in principle the IGM might have ceased to be 
mainly neutral at much higher redshifts, and still produce 
detectable Gunn-Peterson absorption at 2 ^ 6. Note that 
the only way in which we are able to produce a “double¬ 
reionization” is by assuming that stars forming via H 2 cool- 
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Figure 11. Left-hand panel: The Hll filling factor as a function of redshift for models H 2 WeakFB-VMS (solid lines) and H 2 WeakFB 
(dashed lines) in our three cosmologies (see key in right-hand panel). The dotted line shows, for the WMAP PL cosmology only, the 
contribution to the Hll filling factor from stars forming in haloes with virial temperatures below 10‘^K. Escape fractions have been 
chosen so as to produce double reionizations in each model, apart from the Baugh et al. (2005) model (for which we do not compute the 
H 2 WeakFB-VMS model) in which case we simply choose /esc = 1- The escape fractions required are: /esc = 0.25, 0.15 and 1.0 for the 
WMAP RSI, WMAP PL and Baugh et al. (2005) cosmologies respectively. Right-hand panel: The electron scattering optical depth as a 
function of redshift for the same models. 


ing (PopIII stars) are very massive objects producing around 
20 times more ionizing photons per unit mass than for 
a standard IMF. Nevertheless, under such an assumption, 
double-reionization is a natural consequence of our model. 

To allow a rough comparison of our models with obser¬ 
vations of Gunn-Peterson absorption, we have calculated the 
average neutral hydrogen fraction in the IGM, /hi (defined 
as the ratio of volume averaged values of the abundances 
of neutral hydrogen tihi and total hydrogen jih, which is 
the quantity which is most directly comparable to measure¬ 
ments derived from the Gunn-Peterson effect^'^), and used 
this to predict the mean Gunn-Peterson optical depth (top), 
both as functions of redshift. Results for these quantities for 
several of our models (using /esc = 1) are plotted in Fig. 12. 
We also give in Table 3 for each model the (lowest) redshift 
at which the mean neutral fraction drops to /hi = 10“^ and 
10“® respectively, and the redshift for which (tgp) = 3. 

The results for /hi and (tgp) for the different models 
favour the models in the WMAP RSI cosmology if /esc = 1, 
or alternatively /esc appreciably smaller than unity if the 
WMAP PL cosmology is assumed. Note, however, that none 
of our models produce a neutral hydrogen fraction at z = 
6 which is sufficiently high to explain the observed Gunn- 
Peterson optical depths (e.g. White et al. 2003; Fan et al. 
2006) when /esc = 1- (The H 2 model in the WMAP RSI 
cosmology comes closest, achieving tgp = 3 by 2 : = 6.5, but 
by 2 = 6 the optical depth has already dropped to 0.35.) 
The neutral fraction falls very rapidly during reionization 
as can be seen from the fifth and sixth columns of Table 3, 

That is, /hi = J n^iidV/ J undF. This is equivalent to a mass 
weighted average of the local neutral fraction xhi- 


which list the redshifts at which neutral fractions of 10“^ 
and 10~® are reached. Neutral fractions comparable to those 
required by the observational data could be obtained if /esc 
were reduced, at the expense of a slightly lower reionization 
redshift. 

We find that the Baugh et al. (2005) model achieves 
reionization at 2 « 12 for an escape fraction of 100%, while 
adding H 2 cooling to this model raises the reionization red¬ 
shift to 2 « 12.5. The corresponding values of the electron 
scattering optical depth are r = 0.108 and 0.118 respec¬ 
tively. This model is therefore consistent with the WMAP 
measurement of the optical depth, given the rather large er¬ 
rors on that datum. However, the neutral fraction in this 
model is too low at 2 « 6 to explain observations of the 
Gunn-Peterson trough (the optical depth at 2 = 6 for the 
Baugh et al. (2005) model lies in the range 0.01 to 0.1 for 
the various models considered here). While this could be 
rectified by lowering the escape fraction in order to reduce 
the reionization redshift to be closer to 2 « 6, this would 
also reduce the optical depth. 


4.2 CMB Anisotropies 

We have computed the temperature and polarization power 
spectra predicted by our model reionization histories. Fig¬ 
ure 13 shows results for the WMAP PL and WMAP RSI 
cosmologies^^. To calculate CMB temperature spectra, we 
employed the Boltzmann code of Sugiyama (1995), which 

The Baugh et al. ’05 cosmology was not designed to accurately 
match available CMB data and so is not considered here. 
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Table 3. The epochs of reionization for all models considered in this work. Columns 1 through 3 specify input parameters of the models: 
column 1 specifies the cosmological model while column 2 specifies the galaxy formation model used and column 3 specifies the value of 
/esc used. The remaining columns list output quantities. Column 4 lists the redshift of reionization—defined as the redshift at which the 
filling factor in the ionized sphere reionization calculation reaches unity. The first value gives ^reion ^ found directly from our calculations, 
while the second value accounts for the finite upper mass limit of the haloes which we simulate. The true ^reion should lie between these 
two values. The next column gives the electron scattering optical depth (assuming that helium is singly ionized at the same time as 
hydrogen is ionized, and that helium is doubly ionized below 2 = 3). The two values again refer to the calculation without and with 
the correction for the finite upper mass limit of simulated haloes. The next two columns list the redshifts at which the neutral hydrogen 
fraction in the IGM drops to 10“^ and 10~^ respectively, as determined using the IGM evolution model of Benson et al. (2002a). The 
final column lists the redshift at which the mean Gunn-Peterson optical depth reaches a value of 3 


Cosmology 

Model 

/esc 

■2^reion 

T 

%hi=io-2 

^fHI = 10-5 

■^''■GP=3 

WMAP RSI 

Standard 

1.00 

6.9-6.9 

0.061-0.061 

7.36 

7.25 

7.10 

WMAP RSI 

MetFree 

1.00 

7.3-7.4 

0.064-0.068 

7.60 

7.50 

7.38 

WMAP RSI 

TopHeavy 

1.00 

8.7-8.7 

0.079-0.080 

9.02 

8.94 

8.79 

WMAP RSI 

WeakFB 

1.00 

8.1-8.1 

0.073-0.079 

8.50 

8.45 

8.33 

WMAP RSI 

Ha 

1.00 

6.1-6.6 

0.055-0.063 

6.91 

6.69 

6.51 

WMAP RSI 

HaWeakFB 

1.00 

7.8-8.1 

0.074-0.075 

8.44 

8.38 

8.22 

WMAP RSI 

Ha WeakFB 

0.25 

4.4-4.4 

0.048-0.055 

7.60 

7.34 

7.17 

WMAP RSI 

HaWeakFBVMS 

1.00 

10.6-10.6 

0.102-0.102 

8.46 

8.39 

8.21 

WMAP RSI 

HaWeakFBVMS 

0.25 

4.4-4.4 

0.076-0.080 

7.55 

7.36 

7.13 

WMAP PL 

Standard 

1.00 

10.8-10.8 

0.141-0.141 

11.38 

10.74 

10.54 

WMAP PL 

MetFree 

1.00 

11.4-11.4 

0.157-0.157 

11.92 

11.20 

10.93 

WMAP PL 

TopHeavy 

1.00 

15.3-15.3 

0.213-0.213 

14.54 

14.46 

14.33 

WMAP PL 

WeakFB 

1.00 

14.1-14.1 

0.191-0.191 

13.73 

12.77 

12.42 

WMAP PL 

Ha 

1.00 

9.2-10.8 

0.131-0.141 

7.87 

7.79 

7.73 

WMAP PL 

Ha WeakFB 

1.00 

14.1-14.1 

0.206-0.206 

13.87 

12.56 

12.31 

WMAP PL 

Ha WeakFB 

0.15 

6.1-8.3 

0.095-0.107 

7.99 

7.85 

7.84 

WMAP PL 

HaWeakFBVMS 

1.00 

22.2-22.2 

0.299-0.299 

13.77 

12.71 

12.40 

WMAP PL 

HaWeakFBVMS 

0.15 

6.9-8.8 

0.183-0.196 

7.98 

7.85 

7.84 

Baugh et al. (2005) 

Standard 

1.00 

11.9-11.9 

0.108-0.108 

7.97 

7.44 

7.39 

Baugh et al. (2005) 

Ha 

1.00 

12.5-12.5 

0.118-0.118 

7.55 

7.23 

6.87 

Baugh et al. (2005) 

HaWeakFB 

1.00 

14.4-14.4 

0.144-0.144 

8.62 

8.19 

7.76 




Figure 12. The mean mass-weighted Hi fraction (left-hand panel) and Gunn-Peterson optical depth (right-hand panel) as functions of 
redshift for various models in our three cosmologies (see key in figure for colour and line type coding). The model results plotted here all 
assume /esc = 1- The dashed black line in the right-hand panel indicates the fit to the mean optical depth reported by Fan et al. (2006). 
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has been used to calibrate the accuracy of publicly avail¬ 
able codes, e.g., CMBFast (Seljak et al. 2003). The mod¬ 
els shown in the figure have as — 0.9 for the WMAP PL 
cosmology and cts = 0.83 for the WMAP RSI cosmology. 
Furthermore, all models are shown for /esc = 1, except for 
the model labelled “grad:/esc = 0.15” for which we plot the 
H 2 WeakFB —VMS model with /esc = 0.15 (which exhibits a 
double reionization). For comparison, we also show the CMB 
anisotropies predicted for a toy model of instant reioniza¬ 
tion, in which the IGM instantaneously becomes fully reion¬ 
ized at Zreion = 17, choseu to give the optical depth r = 0.17 
that best fits the WMAP data. 

For the WMAP PL cosmology, it is clear that the 
H 2 WeakFB — VMS model with /esc = 1.0 is ruled out from 
the TT-spectrum. The reason for the failure of this model 
is the combination of the large optical depth (which sup¬ 
presses CMB temperature anisotropies on small scales) and 
as ~ 0.9. One could try to make this model agree with the 
TT data by increasing as by ~ 20%. However, this increase 
leads to galaxy formation at higher redshifts and hence a 
larger optical depth. Thus, while increasing as helps to rec¬ 
oncile this model with the high-£ part of the TT spectrum, 
it produces a very large polarization signal at small i, which 
is inconsistent with that measured by WMAP. 

The H 2 and Standard models in the WMAP PL cosmol¬ 
ogy also do not agree well with the WMAP TT-spectrum, 
due to their optical depths being too low. 

It is very interesting that the gradual reionization model 
(the H 2 WeakFB — VMS model with /esc = 0.15; light 
blue line in Fig. 13) shows a unique behaviour in the EE- 
spectrum compared with the H 2 WeakFB and instant reion¬ 
ization models, while the optical depths of these three mod¬ 
els are similar. In particular, this model shows a feature (a 
bump) aX I ~ 20. The location of the £ « 8 peak and the 
width of this peak are slightly different from other models of 
comparable optical depth (i.e. the width is slightly greater 
and the peak location is shifted to the right). These differ¬ 
ences are due to the fact that the reionization process occurs 
over a longer period in this model. On the other hand, there 
are only small differences in TT and TE spectra. We ex¬ 
pect the EE-spectrum to be the best probe of the details of 
reionization such as the duration and the epoch (Holder et 
al. 2003). 

None of the models in the WMAP RSI cosmology fits 
the WMAP CMB data, since we employ the as normaliza¬ 
tion here (i.e. as = 0.83 in all of the RSI models plotted), 
and the optical depths of all of the models are very small. 

The calculations here do not include the effects of patch¬ 
iness in the reionization process. Patchy reionization can 
produce secondary anisotropies in the CMB on small angu¬ 
lar scales, I > 3000 (Paper I). These anisotropies must be 
studied in detail using numerical simulations which provide 
the spatial distribution and peculiar motions of the ionized 
gas. We defer this work to a future paper. 


5 DISCUSSION 

We have described the most detailed and physically realistic 
semi-analytic calculation of the reionization of the Universe 
that has been carried out to date. These calculations rely on 
a physical model of galaxy formation which includes cooling 


through Compton and molecular hydrogen channels, a phys¬ 
ically motivated model of star formation and feedback, and 
the feedback exerted on galaxy formation by the process of 
reionization itself. 

We find that, in the WMAP PL cosmological model, the 
Universe can be reionized by 2 ; « 14--15 either if early gener¬ 
ations of stars are preferentially very massive or if feedback 
in the high-redshift Universe is much weaker than it is to¬ 
day. A direct measurement of the reionization redshift is not 
yet possible. Instead, we can compute, from our reionization 
histories, optical depths for electron scattering which can be 
compared directly with the measurement of this quantity 
from WMAP. Reionization at these redshifts produces op¬ 
tical depths in the range r = 0.19-0.21. This is consistent 
with the WMAP result of r = 0.17 ± 0.04. In the WMAP 
RSI cosmology however, reionization cannot occur before 
2 « 9, due to the much smaller number of massive haloes 
which are able to form at high redshifts in this model. As 
such, the WMAP RSI cosmological model is barely consis¬ 
tent with the WMAP result for r, achieving at most r « 0.1. 
This suggests one of three possibilities: 

(i) the WMAP estimate of r is too high (either due to 
systematic effects or a simple statistical fluctuation); 

(ii) the RSI cosmology is incorrect (perhaps due to prob¬ 
lems in determining the normalization of the Lyman-a forest 
power spectrum; see, for example, Viel, Weller & Haehnelt 
2004) 

(iii) our understanding of galaxy formation, or the CDM 
hierarchical model itself, are incorrect at high redshifts. 

Option (iii) is almost certainly true to some degree— 
there are clearly large gaps remaining in our understand¬ 
ing of galaxy formation. We include within option (iii) the 
possibility that some neglected galaxy formation physics, 
for example a contribution to the ionizing emissivity from 
super-massive black holes forming in high redshift galaxies 
(Ricotti & Ostriker 2004), affects our results. However, in 
the WMAP RSI cosmology, the basic problem is that there 
are simply too few haloes existing at high redshifts to allow 
enough galaxies to form and so cause reionization. Only if 
100% of the baryons which have the potential to cool in high 
redshift haloes are allowed to turn instantly into stars can a 
sufficiently large r be attained in the WMAP RSI cosmology 
(Haiman & Holder 2003; Somerville, Bullock & Livio 2003). 
Within our models, even with rather extreme assumptions 
about the strength of feedback, the IMF of early genera¬ 
tions of stars, and the escape fraction of ionizing photons, 
we are unable to achieve reionization prior to 2 = 9 for the 
WMAP RSI cosmology. Thus, the WMAP RSI cosmology 
seems hard to reconcile with our understanding of galaxy 
formation in the CDM hierarchical model. 

The recent three year data release from the WMAP 
satellite confirms (i) above. Spergel et al. (2006) report r = 
0.093 ± 0.029. This optical depth is well within the reach 
of our models with /esc = 1, although physically motivated 
escape fractions (see Appendix B) would result in models 
still struggling to reach such large optical depths. 

Using the model of Baugh et al. (2005), which was 
designed to match the properties of low and high-redshift 
galaxies, we find that reionization can occur by « « 12 for 
/esc = 1, producing an optical depth of r = 0.11 which is 
marginally consistent with the WMAP data. This optical 
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WMAP PL WMAP RSI 




Figure 13. CMB anisotropy power spectra for the WMAP PL (left-hand panels) and WMAP RSI (right-hand panels) cosmologies. The 
power spectra of temperature fluctuations (TT) are shown in the upper panels, while the middle and lower panels show polarization 
power spectra (EE in the middle panels and TE in the lower panels). The circles with error bars show observational data from the WMAP 
satellite (Hinshaw et al. 2003; Kogut et al. 2003). The dashed black lines (labelled instant) show results for a simple toy model in which 
reionization is instantaneous, and occurs at a redshift selected to best fit the WMAP data. The coloured lines show the predictions for 
various models from this work, assuming /©sc = 100% unless indicated otherwise. The cyan lines (labelled grad) show results for the 
H 2 WeakFB-VMS model with an escape fraction of 15%, which exhibits gradual reionization. 
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depth can be increased to r = 0.12 by including H 2 cooling 
(which does not affect the ability of this model to match local 
galaxy data), and can be increased further if the strength 
of feedback is reduced (which would affect the ability of the 
model to match local galaxy data). 

Recently, the Boomerang experiment has provided 
an improved measurement of the cosmological parameters 
(MacTavish et al. 2005). We have computed reionization 
models using these revised cosmological parameters (specif¬ 
ically, we use the “CMBall-|-B03-l-LSS” parameter set from 
MacTavish et al. (2005), as listed in Table 1). We find that 
at high redshifts (> 20) the Boomerang cosmology produces 
a lower Hll filling factor than the WMAP PL cosmology 
(for all of our different galaxy formation model parameter 
sets). However, this difference becomes less towards lower 
redshifts, and we find that the epoch of reionization in the 
Boomerang cosmology is usually very close to that for the 
corresponding model in the WMAP PL cosmology. The one 
exception is for the H 2 WeakFB-VMS model, in which reion¬ 
ization occurs significantly later in the Boomerang cosmol¬ 
ogy (since reionization happens early in the H 2 WeakFB- 
VMS model, the differences between the Boomerang and 
WMAP cosmologies are more pronounced). Nonetheless, 
double reionization can still occur in our model, even in 
the Boomerang cosmology. We conclude, therefore, that the 
conclusions we present for the WMAP cosmology can also be 
applied to the revised cosmological parameter set suggested 
by the Boomerang data. 

Throughout this work we have used the values of as 
which best fit the WMAP data {as = 0.90 and 0.83 for 
the PL and RSI models respectively, Spergel et al. 2003). 
However, the value of as remains somewhat uncertain, with 
plausible values lying between about 0.7 to 1.1. In particu¬ 
lar, Spergel et al. (2006) report as = 0.76 ± 0.05. We have 
re-run our WMAP PL models using these two extreme val¬ 
ues of as- Results for the ionized filling factor are shown 
in Fig. 14. Not surprisingly, lowering crs reduces the red- 
shift of reionization, while increasing as increases it. The 
effect is greater for models in which reionization happens 
at higher redshifts, since here those haloes which are able 
to form galaxies and emit ionizing photons correspond to 
very rare fluctuations in the initial Gaussian density field. 
As such, their abundance is very sensitive to the normaliza¬ 
tion of that density field. Typically, models which reionize 
at 2 « 10 have their epoch of reionization reduced/increased 
by Az « 1 for a reduction/increase of 0.2 in as, while mod¬ 
els with reionization at 2 « 20 experience a Az of 4-5 for 
the same change in as (these values of Az do depend also on 
the parameters of the galaxy formation calculation). Thus, 
as can have a large effect on predictions of the epoch of 
reionization. 

These results show that, even for the low as reported 
by Spergel et al. (2006) sufficiently early reionization can 
be achieved (i.e. to match the optical depth also reported 
by Spergel et al. 2006) providing /esc = 100%. 

Most of the results we present here are based on the 
assumption that the fraction of ionizing photons which es¬ 
cape from galaxies into the IGM is /esc = 1. In the local 
Universe, and at 2 ~ 3, observational estimates give much 
smaller values. If /esc were low at higher redshifts also, then 
our estimates of 2reion would be substantially too high. A 
value of /esc =0.1 would reduce the reionization redshift in 



Figure 14. The Hll filling factor as a function of redshift in 
various models (distinguished by line colour—see key in figure) 
for different values of erg (distinguished by line type—see key in 
figure). We show results for the WMAP PL cosmology only— 
results for other cosmologies are qualitatively similar. All results 
are shown for an escape fraction of 100%. 


even our best models to 2 < 10. In Appendix B we present 
reionization histories computed using a simple physically- 
motivated model for the escape fractions. Clearly, a better 
theoretical understanding of /esc is crucial if we are to use 
measures of the reionization epoch to put interesting con¬ 
straints on the process of galaxy formation. 

We have found that extended periods of partial reioniza¬ 
tion and double reionization can occur in models in which 
early generations of stars formed via H 2 cooling are very 
massive (with enhanced ionizing luminosities relative to a 
standard IMF). Such models are able to delay full reioniza¬ 
tion until 2 « 6-7, while achieving an electron scattering 
optical depth large enough to be consistent with the results 
from WMAP. While a more extensive study of such mod¬ 
els would be required to determine whether or not they can 
match all of the current observational constraints, they ap¬ 
pear to be a promising route to reconciling the apparently 
discrepant estimates of the reionization redshift from the 
GMB polarization measured by WMAP and from the Gunn- 
Peterson absorption detected in the spectra of 2 ~ 6 quasars. 
It should be noted, however, that these models rely on the 
assumption of very weak feedback from supernovae. Within 
GDM cosmogonies, it has long been argued that feedback 
from supernovae must be very effective in order to reduce 
the number of low-luminosity galaxies at low redshifts which 
are otherwise produced (see Benson et al. 2003a and refer¬ 
ences therein). Thus, if double reionization (or an extended 
period of partial reionization) is required to meet observa¬ 
tional constraints, then it appears that the strength of su¬ 
pernova feedback must vary with redshift, or else some other 
feedback mechanism must act at low redshifts to reduce the 
number of faint galaxies formed. This underlies the urgent 
need for the inclusion of physical recipes for feedback into 











The Epoch of Reionization 21 


semi-analytic models of galaxy formation, rather than the 
phenomenological feedback rules currently employed. 

In conclusion, the WMAP RSI cosmology is very difB- 
cnlt to reconcile with our current understanding of galaxy 
formation (although a wide parameter space of models re¬ 
mains to be searched) if the optical depth were as large as 
suggested by the first year WMAP data. The WMAP PL 
model, in contrast, is able to achieve a sufficiently high opti¬ 
cal depth with reasonable assumptions about the nature of 
galaxy formation at high redshifts. The three year WMAP 
data imply an optical depth significantly lower than that 
found in the first year data (i.e. r = 0.093 ± 0.029; Spergel 
et al. 2006). Reionization can then be achieved sufficiently 
early in both the WMAP PL and WMAP RSI cosmologies, 
providing /esc is close to 1 and ag is not too low. We will 
present a more detailed analysis of the implications of the 
3-year WMAP data in a follow-up paper. 

In principle, the value of the optical depth to reion¬ 
ization is one of the few strong constraints on galaxy for¬ 
mation at high redshifts. The uncertainty mentioned in op¬ 
tion (i) has already been reduced by the WMAP three year 
data.Thus, the ability to reionize the Universe sufficiently 
early will now become one of the key tests that any model 
of galaxy formation must pass before it can be considered 
to be viable. 
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APPENDIX A: REVISED COOLING MODEL 

In this Appendix, we discuss how we incorporate cooling by 
molecular hydrogen and Compton cooling into our model 
and how we apply the effects of the filtering mass. 


Al Cooling by Molecular Hydrogen 

We incorporate cooling due to molecular hydrogen in our 
model following the algorithms of Yoshida et al. (2003a). 
Below, we give a brief account of our calculation of the H 2 
contribution to the cooling rate. The reader is referred to 
Yoshida et al. (2003a) and references therein for a full de¬ 
scription. 

We begin by estimating the fraction, /h 2 , of molecular 
hydrogen which would be present if there were no back¬ 
ground of H 2 -dissociating radiation from stars. As shown by 
Tegmark et al. (1997), at the redshifts of interest here, H 2 is 
formed mostly via the formation of H” through the reaction 
H J- e“ —> H”, with rate coefficient 

fcn- = 1.83 X 10~“(T/K)'’ ®® cm^s"^ (Al) 

(Hutchins 1976). The rate of formation of H” depends on 
the abundance of free electrons, which depends on recombi¬ 
nations -be” ^ H + 7 with rate coefficient 

fci = 1.88 X 10 "“(r/K)"° ®^ cm®s"^ (A 2 ) 

(Hutchins 1976). (We ignore photoionization and collisional 
ionization of hydrogen in this part of the calculation, since 
cooling by H 2 is only important in low-temperature haloes 
at high redshifts.) Tegmark et al. (1997) solved the rate 
equations for electron recombination and H 2 formation to 
obtain a characteristic H 2 abundance /h 2 ,c for gas with hy¬ 
drogen number density nn and temperature T at redshift 
a: 

/ff 2 .c = 3.5 X lO-^Tgi “ [1 + ( 7.4 X 10 ®/nH,i) 


x(l-b 2 )^'^® exp(—3173/[1-b a])] (A3) 

where nn,! is the hydrogen number density in units of cm“® 
and Ts is temperature in units of lOOOK. (The above result 
depends on the values of various rate coefficients, and is con¬ 
sistent with the rate coefficients given in eqns. Al, A2 and 
A 8 .) The factor in square brackets corresponds to photodis¬ 
sociation of H~ by CMB photons, but is important only for 
2 > 100. For our own estimate of the H 2 abundance in the 
absence of a photodissociating background from stars, we 
then use 


/h2,0 — Uiax [/H2,prim,/ h2,c] , (-^ 4 ) 

where /H 2 ,prim is the primordial H 2 fraction, which we take 
to be 3.5 X 10“’^ (Anninos & Norman 1996) (although our 
results are highly insensitive to this choice of value). 

Next, we account for the effects of a background of H 2 - 
dissociating radiation from stars. If we assume equilibrium 
between H 2 production via H“ and H 2 destruction by pho¬ 
todissociation, we obtain an H 2 abundance 

y"H2,eq — /^diss7 (-^ 3 ) 

We calculate the electron fraction, Xe, appearing in 
eqn. (A5) from 

^e,0 / A 

= -—- 7 - 1 -, (A 6 ) 

1 + 3:e,0^1/2/'^rec 

based on eqn.( 6 ) from (Tegmark et al. 1997), where 
Tree = l/(fcinH). (A7) 


This estimate of Xe assumes that recombinations have been 
proceeding in gas at its present density and temperature for 
a time ti /2 equal to half the age of the Universe at redshift 
2 , starting from the primordial post-recombination electron 
fraction Xe,o- We calculate Xe,o using recfast (Seager, Sas¬ 
selov & Scott 2000). Eqn. (A 6 ) accounts crudely for the age 
of the halo, and therefore the time available for recombina¬ 
tions. The electron fraction is never allowed to fall below 
the value corresponding to collisional ionization equilibrium 
at the specified temperature. Note that, post-reionization, 
the electron fraction should in fact tend towards the equi¬ 
librium value for photoionization equilibrium (as opposed to 
collisional ionization equilibrium). Our current model does 
not account for this fact, but this is unimportant here as the 
determination of electron fractions in the post-reionization 
Universe does not affect our determination of the epoch of 
reionization. 

For the photodissociation rate in eqn. (A5), we assume 
fcdiss = 1.1 X 10 ^^J 2 lE 5 hield S (A 8 ) 


(Machacek, Bryan & Abel 2001), where J 21 is average flux 
in units of 10“^^ergs/cm^/s/Hz in the Lyman-Werner bands 
(with photon energy around 12.87 eV), and is computed 
from the radiation background tracked by our model of IGM 
evolution (Benson et al. 2002a). We account for the effects of 
self-shielding using the algorithm of Yoshida et al. (2003a). 
Specifically, we estimate the shielding factor 


Eshieid = min 


A^H2 

10i^cm~2 J 


(A9) 


which depends on the H 2 column density Ah 2 - In order to 
avoid the computational cost and complexity of solving it¬ 
eratively for the radial abundance profile of H 2 within each 



halo, we estimate this shielding factor using the abundance 
/h 2 ,o computed in the absence of a photodissociating back¬ 
ground, using 

A^H2 = C'/h2,0 WH.vir -Rvir (AlO) 

with C = 0.2, where i?vir is the virial radius of the dark 
matter halo and nn.vir is the mean hydrogen density within 
7?vir. (For a more detailed justification of this approach, see 
Yoshida et al. (2003a).) 

We then calculate the final H 2 abundance used in our 
cooling calculation from 

/h 2 = min[/H2,o,/H2,eq]exp (^20k ) ’ 

We have introduced an exponential cut-off in the molecu¬ 
lar hydrogen fraction at temperatures corresponding to the 
dissociation energy of H 2 , to account for collisional dissoci¬ 
ation of H 2 . (Tegmark et al. (1997) ignored this because 
they were interested only in the relatively low temperature 
regime. The exact form and position of the cut-off are unim¬ 
portant, as the cooling function is dominated by atomic pro¬ 
cesses at the temperatures for which collisional dissociation 
is important.) 

Finally, we compute the cooling timescale due to molec¬ 
ular hydrogen using^® (Tegmark et al. 1997) 

rH2 = 48,200 yr(^l + ^j 

xexp(512K/^)(/H2nH,l)■^ (A12) 

The net radiative cooling timescale for gas, TcooI is then given 
by 

Yooi = (A13) 

where Tother is the cooling timescale due to all other cooling 
mechanisms (e.g. atomic processes, Compton cooling). 

Yoshida et al. (2003a) recommend the inclusion of a 
“dynamical heating” term which accounts for the accretion 
of new material onto a halo. This accretion is assumed to 
heat the gas already in the halo (as the kinetic energy of 
the infalling material is converted into thermal energy), re¬ 
ducing the rate at which it can cool. Yoshida et al. (2003a) 
find that the inclusion of this term in their analytic model 
of gas cooling in haloes produces good agreement with the 
results of their gas-dynamical calculations. We have experi¬ 
mented with including such a term in our calculations. We 
note, however, that our semi-analytic model already effec¬ 
tively includes such dynamical heating. Whenever a halo’s 
mass doubles we consider the halo to have reformed (Cole 
et al. 2000). The temperature of the gas in the halo is then 
recomputed, and the cooling clock is reset to zero. It is sim¬ 
ple to show that this results in an effective heating rate for 
the gas which is comparable to that supplied by the dynam¬ 
ical heating term of Yoshida et al. (2003a). For example, 
following the approach of Yoshida et al. (2003a) and taking 
the virial temperature of a halo to be given by 

T = 7 M^/®(l-b 2 ), (A14) 


Using the more recent estimate of the H 2 cooling function from 
Galli Palla (1998) makes only small differences to our results. 
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where 7 is a normalization constant, we find by differentia¬ 
tion that 



(1 + ^) 
Mi/3 


^ M +-M^— 
K 2 1 + z 


(A15) 


(Note that Yoshida et al. (2003a) include only the first 
of the two terms on the right-hand side of this equation. 
However, we expect them to be of similar magnitude.) In 
the semi-analytic calculation, the virial temperature is also 
given by eqn. (A14) and is updated after each halo reforma¬ 
tion event. Consequently, when averaged over a halo lifetime, 
the amount of dynamical heating in the semi-analytic model 
must approximately equal that proposed by Yoshida et al. 
(2003a). 

Not surprisingly therefore, we find that the explicit in¬ 
clusion or otherwise of a dynamical heating term (of the 
form proposed by Yoshida et al. 2003a) makes no signifi¬ 
cant difference to the results of our calculations. 

While in primordial plasma H 2 must form in the gas 
phase via the H“ channel described above, in enriched gas 
there is the possibility of H 2 formation on the surfaces of 
dust grains. While we do not include this formation mech¬ 
anism in our cooling model we have checked that it makes 
no significant difference to the evolution of the filling factor. 
To do this, we used the methodology of Cazaux & Spaans 
(2004) to compute the ratio of H 2 formation rates on dust 
grains and in the gas phase. We then scaled the abundance 
of H 2 in our cooling model appropriately and computed a 
selection of our models using this revised cooling model. We 
found no significant change in the evolution of the H 2 filling 
factors. 


A2 Compton Cooling 

The Compton cooling timescale does not explicitly depend 
on the gas density, but does depend on the electron abun¬ 
dance, and therefore on the ionization state of the gas. The 
electron abundance, in turn, does depend on the gas density, 
due to the density dependence of the electron recombination 
timescale. 

Prior to halo formation, the electron abundance in the 
IGM gas is computed from our model of the ionization state 
of the IGM. Prior to IGM reionization, the electron abun¬ 
dance is low, corresponding to that left over from recombina¬ 
tion. Post-reionization, the electron fraction becomes large 
as the Universe is almost entirely ionized. During and af¬ 
ter halo formation, the electron fraction will move towards 
the appropriate equilibrium value as described in §A1. The 
Compton cooling timescale is then computed using the re¬ 
combined electron abundance from eqn.(A 6 ), or the equilib¬ 
rium electron abundance, whichever is the larger. 


A3 Filtering Mass 

The filtering mass is similar to the Jeans mass, but is the 
relevant quantity for describing the effect of pressure on the 
gravitational collapse of gas in an expanding Universe with 
a time-varying temperature. It is defined by 

Mf — (47r/3)po(27ra/feF)^, (^16) 

where po is the mean density of the Universe, a is the ex¬ 
pansion factor and 
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1 


1 


dt" 


,(A17) 


where D{t) is the linear theory growth factor, H{t) is the 
Hubble variable and kj (t) is the Jeans wavenumber, defined 
in the usual way. The collapse of gas into halos with masses 
below Mp is greatly suppressed relative to the case of zero 
pressure. 


Table Bl. Ionizing luminosity-weighted escape fractions for 
galaxies at two different redshifts (z = 5 and 11) are given for 
all cosmologies and models. Redshift is specified in the first col¬ 
umn, with cosmology and model specified in columns 2 and 3. 
Column 4 specifies the escape fractions for Hi ionizing radiation 
when dust is ignored, while column 5 specifies the same escape 
fractions when dust extinction is included. 


APPENDIX B: ESCAPE FRACTIONS 

Here, we apply simple, physical models to predict the escape 
fraction for each galaxy in our model. We follow Paper I and 
use the model of Dove & Shull (1994) (and the develop¬ 
ment thereof carried out in Paper I) to compute the amount 
of ionizing radiation absorbed by the gas in each galaxy. 
The resulting escape fractions are shown in Table Bl and in 
Fig. Bl. 

The mean escape fraction for Hi ionizing photons for 
galaxies due to absorption by neutral gas are found to be a 
few percent. These numbers agree reasonably well with the 
limited observational determinations of escape fractions, e.g. 
Leitherer et al. (1995) find that less than 3% of ionizing pho¬ 
tons escape from low redshift starburst galaxies, and similar 
results have been obtained by Tumlinson et al. (1999), while 
Steidel, Pettini & Adelberger (2001) estimate escape frac¬ 
tions ~ 10% for Lyman-break galaxies at 2 = 3.4. More 
recently, Bergvall et al. (2006) found an escape fraction of 
4-10% for a local starbursting dwarf galaxy. 

It should be noted that these escape fraction calcula¬ 
tions are extremely simplified—they model the interstellar 
medium (ISM) as a smooth distribution of gas (a clumpy 
distribution could enhance the escape fraction as noted by 
Wood & Loeb 2000) and ignore dynamical effects (such 
as expanding superbubbles as considered by Dove, Shull & 
Ferrara (2000) which significantly lower the escaping frac¬ 
tion). More detailed calculations, employing hydrodynam- 
ical simulations, suggest that high redshift dwarf galaxies 
might have much higher escape fractions (Fujita et al. 2003). 
A better theoretical understanding of escape fractions is a 
crucial missing component of accurate calculations of reion¬ 
ization. 

We estimate the effects of dust extinction on the ioniz¬ 
ing radiation escape fractions by following the methods of 
Paper I. Briefly, we extrapolate the Milky Way extinction 
curve to the wavelengths corresponding to the ionization 
threshold of each species, Hi, Hel and Hell, and compute an 
extinction averaged over galaxy inclinations using the calcu¬ 
lations of Ferrara et al. (1999). The results of including dust 
extinction on escape fractions are shown in Table Bl and in 
Fig. Bl. It is seen that dust can cause severe (frequently an 
order of magnitude) attenuation in the escaping fraction, at 
least in this highly simplified calculation of a smooth ISM. 
Clearly, a more detailed analysis of the effects of dust on 
the escaping fraction is urgently needed if we are to gain a 
better understanding of how the IGM was reionized. 

The effects on the ionized filling factor of using different 
assumptions for the escape fraction are shown in Fig. B2, for 
different models and different redshifts. It is interesting to 
note that, in the WMAP cosmologies, only the TopHeavy 
model in the WMAP PL cosmology is able to achieve reion¬ 
ization by 2 = 6, and then only if the effects of dust on the 
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escape fraction are ignored. In the Baugh et al. ’05 the ma¬ 
jority of ionizing photons are produced during bursts of star 
formation. As such, dust (assumed to be located in galaxy 
disks in our simple model) has little effect on the escape 
fraction. It should be noted that our model for the escape 
fraction is likely to be very unreliable in cases where most 
of the ionizing photons are produced in bursts rather than 
in quiescent disks. 


APPENDIX C: TESTS OF THE MODEL 

In this Appendix we present results of various tests of our 
model. We begin by comparing to the results of Paper I 
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Figure Bl. Ionizing luminosity-weighted mean Hi escape fractions for haloes as a function of halo mass at z = 5 (left-hand panel) and 
2 = 10.6 (right-hand panel) for the Standard model in the WMAP RSI cosmology. The points show the median escape fraction for all 
haloes of this mass, while the error bars show the upper and lower quartiles of the distribution. Filled points include only the effects of 
absorption by neutral gas, while open points also include the effects of dust extinction. 


after which we check that our methodology leads to model 
calculations that are suitably converged. 


Cl Comparison with Paper I Results 

We begin by comparing our current results with those of 
Paper I. Our galaxy formation model, GALFORM, has under¬ 
gone significant development over that used in Paper I. As 
such, we expect some differences in the results obtained. To 
assess the degree of change in our results we have performed 
the same calculations as in Paper I. Specifically, we adopted 
the same model parameters as in that paper and computed 
the Hll filling factor as a function of redshift. The models of 
Paper I use cosmological parameters (flo, Aq) = 0.3, 0.7 and 
(T8 = 0.9, while using similar galaxy formation prescriptions 
and parameters as used in this work. The Paper I models are 
similar to our Standard model, with major differences being 
a weaker feedback, Vhot = 150km/s as opposed to 200km/s 
as used here, and a simpler model for galaxy merging than 
used here (Paper I uses the model of Cole et al. (2000) 
whereas the current work uses that of Benson et al. (2004) 
which gives somewhat longer merger timescales). Further¬ 
more, in Paper I a fraction of all star formation was assumed 
to result in the formation of brown dwarfs, such that the ef¬ 
fective ionizing luminosity per unit mass of stars formed was 
reduced by a factor of 0.65. 

Calculations were performed for both the low and high 
baryon fraction models used in Paper I (fib = 0.02 and 0.04 
respectively). 

As can be seen in Fig. Cl, differences do exist between 
the calculations of this work and those of Paper I. The re¬ 
vised modelling tends to produce slightly earlier reionization 
(by Az « 1). These changes are due to several improvements 


and corrections made to the semi-analytic model and our 
calculation of ionized sphere growth made since Paper I. 

C2 Tests of Model Convergence 

We have examined two possible sources of numerical error 
in our computational method. We outline both below. 

C2.1 Convergence of Ionizing Emissivity 

As noted in §2 we use an iterative procedure to produce a 
self-consistent model of galaxy formation and IGM evolu¬ 
tion. To be self-consistent, the ionizing emissivity as a func¬ 
tion of redshift produced by the galaxy formation model 
should be identical to that used to evolve the IGM in the pre¬ 
vious iteration. We have checked, for all models used in this 
work, that such agreement has in fact been met to within 
the accuracy of the model calculation (which is limited by 
the number of merger trees we are able to simulate). 

C2.2 Upper Limit to Halo Masses 

With finite computing resources, we can only simulate haloes 
with masses up to some maximum value, if we always resolve 
the formation histories of haloes down to the same minimum 
progenitor mass. Our approach is to ensure that the resolu¬ 
tion of the merger trees used in the semi-analytic model is 
sufficient to resolve all haloes in which gas could cool (i.e. 
haloes with virial velocities above 13km/s for models with 
no H 2 cooling or 1.3km/s for models with H 2 cooling). In 
our standard runs, the most massive haloes we consider are 
a factor 1.3 x 10® more massive than the lowest mass haloes. 
The contribution to the ionizing emissivity from galaxies 
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Figure Cl. Comparison of Hll filling factors using the parameters of Paper I. Dashed lines show the results from Paper I (for an escape 
fraction of 100% and no dust absorption). Solid black lines show the results from this work using ionized spheres. Left panel: Model with 
Ob = 0.02. Right panel: Model with Ot, = 0.04. 


in more massive haloes is therefore lost from our calcula¬ 
tion. Figure C2 shows the median Hi ionizing luminosity 
from haloes as a function of mass weighted by the halo mass 
function, dn/d In Mhaio. This therefore gives an indication 
of the relative contribution to the total ionizing luminosity 
from haloes of different masses. At high redshifts, the contri¬ 
bution to the total luminosity declines with increasing halo 
mass (as the halo mass function is steeply declining), so we 
may expect the contribution from higher mass haloes to be 
negligible. At lower redshifts the contribution is an increas¬ 
ing function of halo mass. However, at the highest masses 
shown the relation is quite flat and is expected to turn over 
at still higher masses once the exponential break in the mass 
function is reached. 

We can estimate the effects of the non-inclusion of 
higher mass haloes as follows. We find that, for all models, 
the dependence of the mean ionizing luminosity per halo 
on halo mass can be well represented by a power-law. We 
can use this relation, combined with the known mass func¬ 
tion of dark matter haloes to determine the filling factor 
we would have obtained had we been able to model haloes 
of arbitrarily large mass. We have checked the accuracy of 
this correction by repeating some of our calculations includ¬ 
ing haloes up to a factor of 10 more massive than in our 
standard calculations^®. We find that our correction does 
indeed overestimate the filling factor (as described above), 
but is typically accurate to better than 10%. The correction 


For even larger masses this relation flattens (see, for example, 
Figure 8 of Benson et al. (2000) which shows that, for high mass 
haloes, galaxy formation becomes gradually less efficient as halo 
mass increases), but this will merely make our estimate of the 
effects of ignoring high mass haloes a conservative one. 

While this is possible, the computational time involved makes 
it impractical to do this for all of our models. 
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Figure C2. The Hi ionizing luminosity of haloes of given mass 
weighted by the halo mass function, dn/d In Mjiaio • The points 
show the median value at fixed mass while the errorbars enclose 
80% of the distribution. Results are shown for model H 2 in the 
WMAP RSI cosmology. Filled points correspond to .2 = 6.2 while 
open points correspond to 2 : = 10. The highest halo masses shown 
correspond to the most massive haloes included in our galaxy 
formation calculations at the specified redshift. 


is smaller for higher redshifts (due to the ever more rapidly 
declining halo mass function found at those redshifts), and 
is therefore less important prior to reionization (which is the 
period of interest in this work). 

The relation between halo mass and ionizing luminosity 
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Figure C3. The median number of hydrogen atoms ionized by 
sources in a dark matter halo as a function of halo mass. Results 
are shown for the Standard model in the WMAP PL cosmology. 
Solid lines indicate the median number for redshifts 5.0, 6.2, 8.2, 
10.7, 12.5, 15.4 and 19.8 (increasing from right to left as indicated 
by the arrow). The dashed line shows the relation Nn oc ® for 
comparison. 


which arises from our semi-analytic modelling of galaxy for¬ 
mation is non-trivial as can be seen in Fig. C3. While it is 
well represented by a power-law for the highest mass haloes 
which we consider (A^h oc for haloes in the mass 

range 10^^ to Mq) it steepens for lower mass haloes, 

reaching A^h oc depending on redshift. The form of 

this relation arises from the various physical ingredients (e.g. 
cooling times, supernovae feedback etc.) which go into our 
calculations. Other authors (e.g. Haiman & Holder 2003; 
Wyithe & Loeb 2003; Furlanetto, Zaldarriaga & Hernquist 
2004) have used highly simplified assumptions (such as that 
the number of hydrogen atoms ionized is simply proportional 
to the halo mass; Furlanetto, Zaldarriaga & Hernquist 2004) 
which may be expected to give significantly different results. 

We find that the loss of the highest mass haloes af¬ 
fects our determination of the epoch of reionization signifi¬ 
cantly only for models including H 2 cooling in the WMAP 
PL cosmology. Figure C4 illustrates the effects of missing the 
highest mass haloes in such a model (red lines). For models 
without H 2 cooling our calculations span a sufficient range 
of halo masses to capture nearly all of the ionizing emissiv- 
ity, while for the WMAP RSI cosmology the space density of 
high mass haloes is so small that even models including H 2 
cooling probe a sufficient range of halo masses. The typical 
correction to the Hll filling factor in our models is only 3-5% 
for 2 ; > 10 in models including H 2 cooling. For 2 « 6 the 
correction reaches a factor of 2-3 for models including H 2 
cooling. In models which ignore H 2 cooling the correction is 
less than 1% for 2 > 6. Figure C4 also shows results for the 
model which is minimally affected by the lack of high mass 
halos (Standard model in the WMAP RSI cosmology; blue 
lines) and for a model in which the effect of missing high 



Figure C4. The Hll filling factor as a function of redshift for the 
“H 2 ” model in the WMAP PL cosmology is shown by red lines. 
The solid line shows the result direct from our model calculation 
while the dotted line shows the effects of correcting for the lack of 
the highest mass haloes in our calculation. The horizontal dashed 
line indicates a filling factor of unity. Note that this model has 
the largest fractional correction to the Hll filling factor out of 
those listed in Table 2. For comparison, we show results from the 
Standard model in the WMAP RSI cosmology (blue lines) which 
is least affected by the lack of high mass haloes and from the H 2 
model in the WMAP RSI cosmology (green lines) in which the 
effect of missing the highest mass haloes is more typical. 

mass halos is more typical (model “H 2 in the WMAP RSI 
cosmology; green lines). 

When quoting results for the epoch of reionization we 
list the value computed directly from our calculations along 
with the value obtained after the above correction has been 
applied. The true epoch of reionization should then lie be¬ 
tween these two values. 


APPENDIX D: CONDENSED FRACTIONS 

Figures D1 and D2 show condensed fraction (i.e. the mass 
of material condensed into the galaxy phase divided by 
HbMhaio/Ho) as a function of halo mass Mhaio for redshifts 
2 = 5 and 10 respectively. Columns indicate the cosmology 
(see labels in upper row) while rows indicate the model (see 
labels at right-hand edge of each row). In each panel, the 
vertical dashed magenta line indicates the virial tempera¬ 
ture corresponding to the feedback parameter Vhot while 
the vertical dashed cyan line indicates the filtering mass. 
Red points indicate the condensed fraction (note that the 
vertical scale is logarithmic). 

The effects of including H 2 cooling can be clearly seen 
as the presence of a tail of low /cond for low mass halos. 
The effects of weak feedback are also apparent, resulting in 
a small increase in /cond in low mass halos. 
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Figure Dl. The condensed fraction, /cond? defined as the mass of material in a halo which has condensed into the galaxy phase, divided 
by f^b^halo /^0 is shown by red circles. Results are shown for 2 = 5 and for all three cosmological models (one per column as indicated 
by the labels in the upper row of panels), and for all models computed (one model per row as indicated at the right-hand edge of each 
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Figure D2. As Fig. D1 but for 2 = 10. 
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